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INTRODUCTION 

This  report  documents  the  numerical  analysis  for  the  hydrofoil  design 
procedure  given  recently  by  Parkin  and  Fernandez  (Ref.  1).  It  also 
presents  a computer  program  for  the  implementation  of  this  design  method, 
along  with  outputs  from  sample  runs. 

The  design  method  given  in  Ref.  (1)  uses  a linearized  cavity  flow 
theory  to  determine  the  wetted  surface  shape  from  a prescribed  pressure 
distribution  and  to  evaluate  its  overall  performance  parameters,  while 
controlling  the  cavity  thicknesses  at  two  points  along  the  chord.  These 
thicknesses  are  specified  at  the  trailing  edge  and  at  an  arbitrary  point 
between  the  leading  and  the  trailing  edges.  For  further  details  of  this 
design  method,  the  reader  should  consult  Ref.  (1). 

The  computer  program  is  coded  in  FORTRAN  IV  and  has  been  implemented 
to  run  on  IBM  system  370/168.  The  program  contains  about  1200  statements 
and  consists  of  a MAIN  program  and  34  subprograms.  The  main  program  reads 
input  data  and  calls  the  principal  subroutine  SECT6,  which  performs  the 
necessary  computations  for  the  design  method.  The  input  to  the  program 
is  described  in  the  comment  statements  for  the  MAIN  program.  The  program 
listing  and  the  sample  output  are  given  in  the  Appendix.  The  following 
sections  describe  the  implementation  details  for  the  various  subprograms. 
For  those  subprograms  that  are  included  in  the  program  listing,  but  are 
not  described  in  this  report,  the  reader  is  directed  to  References  (2)  and 
(5). 

SUBROUTINE  SECT 6 

SECT6  is  the  major  subro  'ne  for  the  third  design  method.  This 
subroutine  operates  in  two  moae^,  namely,  the  estimation  mode  (NEST=1) 
and  the  design  mode  (NEST*0) . 


f 
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In  the  estimation  mode  SECT6  performs  the  following: 

(1)  Estimate  the  value  of  CT  at  which  y (1)**0  for  a 
prescribed  m.  This  is  the  maximum  permissible 
value  of  for  the  given  m. 

(2)  Estimate  the  value  of  y for  the  prescribed  CT  and  m. 

la 

For  the  general  case  K?*0,  the  above  estimates  are  obtained  by  solving 
for  the  cavity  length  parameter  'a'  by  using  the  appropriate  equation.  These 
equations  are  evaluated  by  FUNCTION  FA.  For  the  particular  case  K=0,  the 
above  esimates  are  calculated  in  SECT6  itself. 

In  the  design  mode,  SECT6  computes  the  cavity  shape,  y(x),  camber  curve, 
n(x)  and  the  slope  of  the  camber  curve,  ri'(x),  along  with  various  performance 
parameters.  In  this  case,  the  designer  prescribes  y from  the  range  of 
y values  obtained  in  the  estimation  mode.  The  following  operations  are 
performed  by  SECT6  in  the  given  order. 

Case  A.  K^O 
Step  1: 

The  transcendental  equation  f(a)=0  is  solved  for  a In  this  step.  The 
expression  f(a)  is  given  in  the  description  of  FUNCTION  FA  for  various 
cases  that  arise  in  the  estimation  and  design  modes.  The  method  for 
obtaining  'a'  consists  of  two  steps.  The  first  step  uses  a sequential  search 
procedure  to  obtain  a^  and  a^  such  that  ' 
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Thus  the  root  of  the  equation  f(a)-0  Is  located  within  the  range 
< a < a^.  In  the  second  step,  the  value  of  'a'  is  obtained  to  within  a 
prescribed  tolerance  using  the  method  of  false  positions  (Regula  Falsi) . 

At  this  point  control  is  returned  to  the  calling  routine  if  the  program 
were  in  the  estimation  mode.  Otherwise,  SECT6  obtains  the  following 
design  parameters. 


Step  2: 


The  parameter  B is  obtained  from  the  equation 


_ <1'm)CL  (e-.-a5)  m f q-2x)P(x; 

K(l+K)  a£(o-ae)  4Tr(l+K)a£  J 


The  value  of  the  integral  in  the  above  equation  is  obtained  in 

the  subprogram  FUNCTION  FA  and  is  transmitted  through  common  as  B^. 

* 

The  parameters  e and  0 in  the  above  equation  are  defined  by 


{6}  = [i/£+l]1/2  ± [✓T-l] 1/2 


where  £,  the  cavity  length  is  obtained  from 


* _ 

The  following  expansion  formulas  for  e and  o are  found  useful  in  obtaining 
various  limiting  values  as  a-*»  (or,  equivalently,  as  K-K>) 


S * /l  {2  + — ~ + 0(^-)} 

4a  a 


E^u'^  + 0(y 
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Step  3; 

The  parameter  D is  obtained  from  the  relation 

(lm)  CL  <5  m f P(x)dx  B 

° " tt  (1+K)  6-ae  4tt(1+K)  J £-x  “2 

0 


The  integral  in  the  above  equation  is  obtained  in  subprogram 
FUNCTION  FA  and  is  transmitted  through  Common  as  B^. 

Step  4: 

Next,  the  value  of  the  parameter  E is  obtained  from 


3 - 


U-m)  CL  (2a2+3)e  + 


a<S 


m 


2tt(1+K) 


£(6-ae) 


2t7(1+K)£ 


/ (£-jx) /x  P(x)  dx 


The  value  of  the  integral  in  the  above  equation  is  transmitted  as 
from  FUNCTION  FA. 

Step  5; 

The  parameter  A is  obtained  from  the  relation 


-2U-B)  CL  A 
A ° TT(l-fK)  (6-ae) 

Step  6: 

Next,  the  value  of  the  drag  coefficient,  Cp,  is  calculated  from 
CD  = 2it(1+K)  £(aB  + |)2 


a 


7 1 


(4) 


(5) 


(6) 


j 


(7) 


J 
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Step  7: 

The  moment  coefficient  is  obtained  from  the  relation 


CL,  - -(1-m)  C.  i3i-Aa^_a2(6,-2ae)  . m f x p(x)  dx 
^ L 8e/l  J 


(8) 


It  should  be  noted  that  the  factor  (36-4ae)  in  the  above  formula  is 

1/2  ■?/? 
of  the  order  a ' while  the  factor  (S-2ae)  is  of  the  order  (1/a  ' ) . 

Therefore  the  later  has  to  be  computed  in  double  precision  in  order  to 

ensure  accuracy  at  large  values  of  a.  The  integral  in  the  above  formula 

for  is  computed  by  the  function  subprogram  BICM. 

Step  8; 

Now  the  contribution  from  the  nose  singularity  to  the  cavity  thickness 

at  the  trailing  edge,  viz.  y (1),  is  calculated  from 

c 


ycU)  = aE 


(9) 


Step  9: 

The  parameter  V..(0)  is  needed  for  the  off-design  calculations.  V (x) 
* * 

is  the  well-behaved  part  of  the  camber  derivative  function  and  is  defined 
as 


E 

V.  (x)  * n'(x)  + — — - , 0<x<2h 

2/jc 

*=  n'(x)  , x>2h 


(10) 


where  and  are  respectively  the  coefficient  of  the  complementary 
function  and  the  cavity  length  at  the  design  conditions.  V^(0)  is 


obtained  from 
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_ (1-m)  CL  5 

Vl(0)  - a - t(1+k)  (6-ae) 


(ID 


Step  10; 

The  cavity  function  y(x)  is  obtained  for  points  along  the  chord  from 
the  relation 


y(x)  * ax  - 2aZ£  AF^(x,a)  + aB[  A(£-x)  - £ tan  ^ A/  < £-x)  ] 


(|  + D)x  + eA(£-x)  + 2^°^  G(x,a;P) 


(12) 


The  function  F^Cx.a)  is  described  in  Ref.  5 and  is  computed  by 


FUNCTION  F1FUNC. 

The  function  G(x,a;P)  is  given  by 


-1 


G(x,a;P)  = | B + — [A(£-x)  + £ tan  ± /x/(£-x)]  B 


2 3 2 


tan  1 A/(£-x)  B^  + £n  (j^) 


(13) 


x 

-I 


P(t)  £n 


1 + /x(£-t)/t(£-x) 


dt 


The  quantities  B^,  B^  and  B^  represent  integrals  involving  the  function 
P(x).  These  are  obtained  in  FUNCTION  FA  and  are  transmitted  through 
COMMON.  The  integral  in  the  above  equation  can  be  rewritten  as 


x x 

j P(t)  £n  1 + A(£-t)/t(£-x)  | dt  = j P(t)  £n(A(£-x)  + A(£-t))  dt 


-if 


P(t)  £n  t dt  - y-  £n(£-x) 


11 


(14) 
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The  first  integral  on  the  right-hand  side  of  the  above  equation  is 
obtained  from  subprogram  FUNCTION  BTX7.  The  seccnd  integral  is  given 
by  BTG5. 

Step  11: 

The  camber  function  p(x)  is  obtained  from  the  relation 


U(x)  * ax  -(|  + D)x  + aB[£  tan  1 /x/  (£-x)  - /x(£-x)  ] - E /x(£-x) 


-1 


m rx 


2u (1+K) 


{y  B^  - y[v'x(£-x)  + £ tan  ^ /x/ (£-x) ]B^  + tan  ^ /x/ (£-x)  B 


c 1 
+ y-  In  (^)  - j P(t)  £n 

1 - /x(£-t)/t(£-x) 

0 

dt} 


(15) 


The  quantities  B^,  B^  and  B^  are  obtained  in  FUNCTION  FA  and  are  transmitted 
through  COMMON.  The  logarithmic  integral  can  be  rewritten  as 


1 

0 

P(t)  £n 

0 

1 - A(£-t)/t(£-x) 

dt 


CT  in  (£//£^)  + 

L 


i i 

1 

2 


P(t)  £n(t-x)dt  - -r 


P(t)  in  t dt 


- | P(t)  £n(/t(£-x)  + /x(£-t))dt 


(16) 


The  first  integral  on  the  right-hand  side  of  the  above  equation  is  given 
by  FTJX,  the  second  by  BTG5  and  the  third  by  BTX7 . 

Step  12: 

The  derivative  of  the  camber  function,  viz.  r)'(x)  is  calculated  in  this 
step,  u' (x)  is  given  by  the  relation 


4 
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n’(x)  - a - f (1  - 2a/x/(£-x))  - D - § — ~2-- 
Z Z /x(£-x) 
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{2  B3  - 2 B1  + 


(17) 


+ ft  __x_ 

' A ✓ /I  \ ' 


2/x(l-x) 

1 p>v>  | P(t)(/t(&-x)  + A(A-t))/TF  dt} 


2 (A_x)  2A(£-x) 

v/ 

The  quantities  B^,  and  are  obtained  through  COMMON  from  FUNCTION  FA. 

The  Cauchy  principal  value  of  the  singular  integral  is  obtained  from  subprogram 
FUNCTION  CPVI. 

Step  13: 

The  overall  design  parameters  like  C , C^,  y (1),  A,  B,  D etc.  are 

M D c 

printed.  The  cavity  shape,  camber  shape  and  the  camber  derivatives  are  also 
printed  as  functions  of  x. 

Case  B:  K=0 

In  the  estimation  mode,  SECT6  performs  Steps  14  through  17.  In  the 
design  mode.  Steps  16  and  18  through  27  are  performed. 

Step  14: 

In  this  step,  SECT6  computes  the  maximum  value  of  C^  for  the  case  K=0. 


This  maximum  permissible  CL  corresponds  to  yc(l)=0.  The  C^ay  is  given 


by 


11- 
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Thus 


1 

| q(x)dx  - 1 

0 

The  first  integral  in  Eq.  (18)  is  calculated  using  BTG9.  The  second 
Integral  can  be  rewritten  as 


1 

* 

0 


q (x)  Zn  dx 

1-& 


1 


2 


j q(x)£n(l-tVx)dx 
0 


q(x)£n(l-x)dx 

0 


(20) 


The  first  integral  on  the  right-hand  side  of  the  above  equation  is 
evaluated  by  FUNCTION  BTG8  and  the  second  by  FUNCTION  FTJX. 

Step  15: 

This  step  computes  the  value  of  m corresponding  to  yc(l)=0  for  a 

prescribed  CT . This  value  of  m is  given  by 

JL 


l/l  - ln(l+/2)  - uT/Ct 


m 


3 Jl  - ln<W2>  - f aiiM*  + I 


q (x)  £n(— ^)  dx 
1-v^t 


(21) 


The  integrals  in  the  above  equation  are  computed  as  described  in  Step  14. 
Step  16: 


This  step  computes  the  quantities  C^  and  C C^  and  C ^ are  defined 
as 


-12- 
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C.  “ 2C.  [{3/2  - An(l+/?)}  4 + “ => 

1 L L iJT 


- (1+2^T)  /(1+v^T)  x "3/4  + — Jln(xrt1/4  + A +/3T)  ] 


x o 

o 


X 

C2  - ci  + f P(x)  Jln(l  + — : 

0 S* 


i->  dx  - ±-  f 

/“  X J 

*/v  r»  * 


r ^T 

- P(x)  Jin  (1  + — ) dx 
0 0 ^ 


j- 

- 4 - -~)  f P(x)  An  dx  . 

2^r  J 1-^ 


The  quantity  C2  can  be  rewritten  as 


l x 

C2  “ — — J P(x)  £n(l  + /x)dx  - j P(x)  An(y^c  + »/x^)  dx  + 

/xo  0 X°  0 

1 1 

- ~ (1  - — ) I P(x)  Jin  x dx  + 4 - — — ) P(x)  Jln(l-x)dx  . (23) 

1 Xo  J L l/x~ 

0 0 o 0 

The  first  and  the  second  integrals  in  the  above  equation  are  computed 
by  SUBROUTINE  BTG8,  the  second  by  SUBROUTINE  BTG5  and  the  third  by  SUBROUTINE 
FTJX. 

Step  17: 

This  step  computes  y for  a given  m in  the  estimation  mode,  y is  given 


<m  C2  - Cx)  xq  xq  + 


Step  18: 


V 


The  value  of  m is  calculated  in  this  step  for  a given  y in  the  design 


mode,  m is  given  by 


C.  + 2ttT(J-  - i ) 

A o Z 2r/x~ 

o 


Step  19: 
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The  value  of  y£(l)  is  computed  from  the  formula 


(l-m)C 

yc(1)  mi  — or" {3v^  * £n(1  + ^ )} 


- fj  < [ | f P<x)  *,.<— > dx} 

J 1 — */ V 


The  first  integral  in  the  above  equation  is  given  by  BTG9  and  the  second 
is  computed  as  described  in  Step  14. 

Step  20: 

The  value  of  a,  the  design  attack  angle  is  computed  in  this  step  from 
the  relation 


2(l-m)CT  f . . . f . 

+ yc(l)  + |ijjr  { — — dx  + P(x)£n(l  - — )dx} 

J JZ  i JZ 


The  first  integral  in  the  above  equation  is  given  by  BTG9.  The  second 
can  be  rewritten  as 


X X 

P(x)  £n(l  - — ) dx  = P(x)  £n(l-x)dx  - -j  P(x)  £n  x dx 

i fir  • 


- | P(x)  £n(l  + >/x)dx 


The  first  integral  on  the  right-hand  side  of  the  above  equation  is  given 
by  FTJX,  the  second  by  BTG5  and  the  third  by  BTG8. 


-14- 
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Step  21; 


The  drag  coefficient  CD  is  computed  in  this  step  from  the  relation 


x 

CD  - ^ UU-m)^  + f | 


P(x)dx  . _ ,,,.,2 

-~^-t — + Try  (1)} 

Jk  c 


(29) 


The  integral  in  the  above  equation  is  given  by  BTG9. 
Step  22: 

The  moment  coefficient  is  obtained  from 

1 


*=  -(1-m)  CL  - m | P(x)  x dx 


(30) 


The  integral  in  the  above  equation  is  computed  by  FUNCTION  BICM. 
Step  23; 


The  value  of  V^(0)  needed  by  the  off-design  calculations  is  computed 


in  this  step.  V^(0)  is  given  by 


V1(0)  = a - 2 (1-m)  -f  . 


(31) 


Step  24: 

This  step  computes  the  cavity  shape  y(x)  for  points  along  the  chord 
for  the  case  K=0.  This  is  given  by 


(1  — m)C_  a It  . tl 

y(x)  - ax  + — {(l+2*^c)x  A+/x  - 2x  - £n(x  + /l+/x) } 


+ yc(l)  Sk  + ^ {Jk  | P^d-t~  - 


P(t)  JLn(l  + — ) dt} 

J t 


1 


(32) 
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The  first  Integral  is  given  by  BTG9  and  the  second  can  be  rewritten 


as 


1 

r 


J P(t)  £n(l  + — )dt  - | P(t)  in  (,71  + v£)dt  - j P(t)in  t dt 

0 * i 0 


(33) 


The  first  integral  on  the  right-hand  side  of  the  above  equation  is  given 
by  BTG8  and  the  second  by  BTG5. 

Step  25: 

In  this  step,  the  camber  shape  of  the  hydrofoil  is  obtained  from 
the  result 


n(x)  * (a  - 


2(1-m)CL 

— )x  - yc(l) 


P(t)dt 

7t 


-fcjp( 


(t)in  ( 


i-£ 

71 


)dt 


(34) 


The  first  integral  in  the  above  formula  is  given  by  BTG9  and  the 
second  can  be  rewritten  as 


P(t)in  ( 


a 


)dt  = | P(t)in( | t-x | )dx  - | 

0 0 
1 

- | P(t)in  (7t  + t/x)dt  . 


P(t)in  t dt 


(35) 


The  first  integral  on  the  right  is  given  by  FTJX,  the  second  by  BTG5 
and  the  third  by  BTG8. 
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r 


n’(x) 


. 2(1-3)  c 

it  CL 


yca> 


+ ~ — {P.v.  f ISmi.  [L(t)dt} 

4ir«^  A-/Z  ^ /t 


(36) 


The  Cauchy  Principal  value  Integral  in  the  above  equation  can  be 
rewritten  as 


’ 


PQQdt 


1 

f P(t)  (/t-JVx) 


t-x 


0 


dt 


(37) 


The  integral  on  the  right  is  evaluated  by  function  CPVI. 

Step  27: 

The  overall  design  parameters  like  K,  C^,  T,  CM,  CD,  etc.  are  outputted 
in  this  step.  The  cavity  shape,  camber  shape  and  the  camber  slope  are  also 
outputted. 

SUBROUTINE  EST 

This  subroutine  computes  the  limiting  values  or  "estimates"  the  parameter 
U for  prescribed  CT  and  other  input  variables.  These  estimates  are  obtained 

1j 

in  two  steps. 

Step  1: 

The  maximum  values  of  CT  are  obtained  at  m=0  and  m=l.  These  are 

Lj 

designated  as  CLq  and  Cj^  respectively.  The  maximum  C^  is  obtained  by 
setting  m to  the  required  value.  The  program  switches  NEST  and  NCLM  are 
set  to  1 and  the  subroutine  SECT6  is  called.  The  maximum  C^  is  obtained 
through  COMMON  upon  return  from  SECT6. 

Step  2: 

The  limiting  value  of  p is  obtained  for  a specified  m and  in  this 
step.  The  switch  NEST  is  set  to  1 and  the  switches  NEM  and  NCLM  are  set 


.. 


. •'  T 
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Co  0.  Various  cases  arise  depending  on  the  value  of  CT . 
Case  a: 


CL  < min  (Cl0  , Clx) 


The  limiting  y values  are  those  corresponding  to  m=0  and  m=l.  These 
are  obtained  by  calling  SUBROUTINE  SECT6. 

Case  b: 


CL  < max  (CLq  , CLl) 


In  this  case  the  diagnostic  message  "Design  CL  is  greater  than  maximum 
permissible  C^"  is  printed  and  the  control  is  returned  to  the  calling 
routine. 

Case  c: 


CLl  £ cL  1 cLo  . 


In  this  case  the  switch  NEM  is  set  to  1 and  the  value  of  m'  in  Eq.  (39) 
is  obtained  from  SECT6.  The  limiting  values  of  p correspond  to  m=m'  and 


Case  d: 


Cl0  * CL  * Cl1 


The  value  of  m'  is  obtained  from  Eq.  (39)  by  setting  NEM«1  and 
calling  SECT6.  The  limiting  values  of  u correspond  to  m=m'  and  m=l. 
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SUBROUTINE  FA 

This  subroutine  computes  the  limiting  values  of  CL  and  y and  also 
calculates  the  appropriate  functions  f(a)  for  the  Iterative  solution  of 
the  equation  f(a)**0. 

Step  1: 

The  parameters  c and  6 are  computed  from 

{*}  - [/I+l]1/2  ± [/I-l]1/2  (3) 

2 

where  £=l+a  . 

Step  2: 

This  step  computes  the  limiting  value  of  C for  a given  m.  This  is 
given  by 


2tt(1+K) 


(1-m)  ~~a  + m / /(f,-x)/x  q (x)  dx 

6-ae  0 


(38) 


where  q(x)  is  the  normalized  pressure  distribution  defined  by 


q(x) 


-POO 


P(x)dx 


Step  3: 


This  step  computes  the  value  of  m corresponding  to  yc(l)=0  for  a 


given  C . This  is  obtained  from 
JL 


Ik.  k 
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02it(1+K)JI/Ct  * T 


A^~x)/x  q (x)  dx  - T 


The  following  integrals  are  evaluated  in  this  step.  The  values  of 
these  integrals  are  made  available  to  other  subprograms  through  COMMON 


computed  by  FUNCTION  BINT1 


B,  is  given  by  FUNCTION  BTG1 


P(x)dx 

<*-x) 


B,  is  obtained  from  FUNCTION  BINT3 


B.  Is  obtained  from  FUNCTION  BINT4 


3.  _ : . ■«*.— 
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B_  - I P(x)  - JIB,  - 2B, 

5 J0  /^o  1 4 


B6  - P(x)dx  « CL 


X 

B,  - f P(x)  In  £5  dx 

7 I,  + /E5 


By  can  be  rewritten  as 


By  « £n(£.)  + | P(x)£n|x-l|dx 

0 

1 

- 2 | P(x)  £n(a>^c  + /£-x)  dx 


The  first  integral  in  the  above  equation  is  c mputed  from  FUNCTION  FTJX, 
and  the  second  by  FUNCTION  BTX7. 


Step  5: 


The  quantities  and  F^  are  computed  in  this  step.  These  are  given 


CT  4a^  i/a 


-1  As 


>n  - (i-x)  Fi(1-a>  + <“  - 11  ta"  <:»  cl 


Jl(6-ae) 


a „ (2a  +3)e  + a6  . a „ 1 r „ _ -1  ,1m  „ 

-ICL  + TB2  - 21  !a  - 1 t,n  (;)}  B5 


\ {a  + n tan-1  <±)}  Bj  + tan_1(i)  B4  - j-  . 
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Fn  “ T ■ 2aB  ■ if (1+K)  "(6-ae)  F1  (1,a) 


L r „ . -1  .I.-,  e+a6 

la  - x tan  (— )}  777 — rr 

if  (1+K)  'a  Jt(o-ae) 


a (2a2+3)e  + a6 


if  (1+K)  SL  6-ae 


The  function  F^(x,a)  is  obtained  from  F1FUNC. 

Step  6; 

The  quantities  Ht  and  F^  are  computed  from  the  equations 

Ht  “ cl  {Fi(1*a)  - “hr — }{a  • /(£_xo)/xo ' 1 tan’1 

o 

. 1 . -1  r—TT-n n r°L  e+a6  ®5i  , r®2  CL  (2a2+3)e  + aS, 

+ — tan  /xo/  («-xo)  }{—  ■)  - 44}  + { 24  " y ‘ Jl(6-ae)  } 

o 

G (x  ' a* P) 

X (a  - A*-  xo)/xQ)  + I {G(l,a;P)  - — - — } . (' 

o 


CT  , 20  r-  F (x  ,a) 

p = t L — 4a  2/a  r^,  q _ _1 2 } 

Ft  tt(1+K)  (6-ae)  1 1U’  ' xo  1 


- {a  - A^-x  )/x  - l tan  * (— ) + — tan  1 A / (Jl-x  ) } — -,1^  . 0 ?ta<~’v 
00  ax  00  2tt(1+K)  x,(6-ae) 

o 


f2a  +3)e  + a6 


2it  (1+K)  £(6-ae) 


Ha  - 


t-x  /x  J 
o o 


The  function  F^(x,a)  is  obtained  from  F1FUNC.  The  function 
G(x,a;P)  is  given  in  Eq.  (13). 
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Seep  7: 


This  step  is  executed  if  the  program  is  in  the  estimation  mode  (NEST=1) , 
This  step  computes  the  appropriate  function  f(a)  used  in  the  iterative 
solution  for  a.  f(a)  is  given  by 


f<*>  * Fn  + 70%  Hn 


The  corresponding  value  of  y is  obtained  from 


* ' r «, + tost  V 


The  control  returns  to  the  calling  routine. 


Step  8: 


This  step  is  executed  if  the  program  is  in  the  design  mode  (NEST=0) . 
The  function  f(a)  is  given  in  this  case  by 


f(a>  - Fn  Ht  - (ft  - f>  Hn 


The  corresponding  value  of  m is  obtained  from 


m = - F tt(1+K)/H 

n n 


The  control  is  now  returned  to  the  calling  routine. 

FUNCTION  CPVI 

The  function  subprogram  CPVI  evaluates  the  following  Cauchy  Principal 
value  integral 
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r 

1 


J m j11  P(t)  g(t,x)dt 
X1 


(55) 


where 


g(t,x)  - (A(i-x)  + ATFO)  A-t  , K^O 


* (/t  + A) 


K=0 


The  value  of  x and  Its  Index  I in  an  array  of  x values  are  inputted 
to  the  subprogram  CPVI.  Various  cases  arise  depending  on  the  value 
of  1. 

Case  A:  1*1 

In  this  case,  we  write 

2 

P(t)  g(t,x)  = a^t-x)  + a2(t,x)  , xx<t<x3 

Since 


P(x)  = 0 for  1=1 


(56) 


The  values  of  a^  and  a2  are  obtained  from  DIFF5.  The  integral  J can  be 
written  as 

Xj+2h 

J = | [a^  + a2(t  - x^) ]dt  + j dt  . (57) 


The  first  integral  is  evaluated  as 
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xx+2h 

[ax  + a2(t-x1)]dt  = 2h(aL  + a2h)  . (58) 

X1 

The  second  Integral  Is  evaluated  by  subprogram  DIFF2. 

Case  B:  1=2 

In  this  case,  we  write. 


P(t)  g (t,x) 


aQ  + a1(t-x) 


+ a2(t-x)‘ 


xx  i t < x3 


(59) 


where  the  coefficients  a^,  a^  and  a2  are  obtained  from  DIFF5.  The  integral 
can  be  written  as 


x,+2h 


"1  • *■“ 

I + >1 + + 1 


(60) 


x^+2h 


The  value  of  the  first  integral  is  given  by  function  F(x2,  x,  x2+2h, 
Sq,  a^,  a2)  where  F is  defined  as 


x -x 


F(xr  X,  xu,  q0,  qv  q2)  = qQ  In  (~^>  + ql(xu~x£) 


q2 

+ 2~[  (xu~x)  ~ (xj-x)  )]  . (61) 


The  second  integral  can  be  evaluated  by  DIFF2. 
Case  C: 


3 i H N-3 


In  this  case,  the  integral  J is  written  as 
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j . f 1 dt  + I"  + «!  + a2(t-x)]dt 

X1  Xi-1 

+ f"  ISiLiii**!  dt  . 

J t-X 


The  first  and  the  third  Integrals  in  the  above  equation  are  evaluated 
by  DIFF2.  The  second  integral  is  given  by  FCx^,  x,  xi+1>  aQ,  a.^  a.^) , 
where  the  constants  a^,  a.^  and  are  defined  by  the  relation 


F(t)  g(t,x)  = aQ  + a^(t-x)  + a2(t-x)  , xi_i^t^xi+i  * (6 


Case  D:  I=N-2 


The  integral  J is  written  as 


J - f 3 F(t)  £'*)dt  + f’1  + »i + V'-”-1' 

*1  *N-3 

> co 

+ J feic  + C1  + c2^t_x^dt 

Vi 

where  (b^,  b^,  b2)  are  obtained  from  the  parabolic  fit 


P(t)  g(t,x)  = bQ  + b^t-x)  + b2(t-x)  , 


£ t £ 


and  (cft,  c, , c_)  are  obtained  from  the  approximation 
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P(t)  g(t,x)  = cQ  + c^t-x)  + c2(t-x)  , xN_2<t<xN 


The  first  integral  in  Eq.  (64)  is  evaluated  by  DIFF2  and  the  second 
and  the  third  integrals  are  given  by  function  F in  Eq.  (61). 

Case  E:  I*N-1 


In  this  case  the  integral  J is  given  by 


J , j”  2 PCtj_g.(.t^x)dt  + | [JL  + bx  + b2(t-x)]dt 

X1  V2 


(66) 


where  (b^,  b^,  b2)  are  obtained  from  the  approximation 


P(t)  g(t,x)  = bQ  + bx(t-x)  + b2(t-x)  , Xjj^^tlXjj 


(67) 


The  first  integral  in  Eq.  (66)  is  given  by  DIFF2  and  the  second  by 
function  F of  Eq.  (61). 

Case  F : I=N 


In  this  case,  the  integral  can  be  written  as 


*-2 


ft 


J = | P-^  + [b1  + b2(t-x)]dt 


(68) 


*N-2 


where  the  constants  b^  and  b2  are  obtained  from  the  approximation 


P(t)  g(t,x)  = b1(t-x)  + b2(t-x)  , x^StSXjj 


(69) 


which  satisfies  the  condition  P(xN)=0. 
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MAIN  LINE  PROGRAM  FOR  THE  THIRD  DESIGN  PROCEDURE  (SECT6 ) 
THIS  PROGRAM  READS  THE  INPUT  DATA  AND  CALLS  THE  PRICIPAL 
ROUTINE  SECT6  FOR  THE  HYDROFOIL  DESIGN  WITH  TIO  POINT 
CAVITY  THICKNESS  CONTROL. 

JUNE  15,  1977 


DICTIONARY  OF  SYMBOLS 
A IS  THE  CAVITY  PARAMETER  LITTLE  A. 

ACAP  IS  THE  'NOSE  SINGULARITY'  PARAMETER  CAPITAL  A. 

A J 1 , A J2 , A J3 , A J 4 ARE  THE  EXPRESSIONS  J 1 , J 2 , J 3 , J 4 W H I C H INVOLVE  THE 
FUNCTIONS  F1(NU,A) ,F2(NU,A) ,L(NU,A) ,AND  T(NU,A)  — SEE  SECT5. 

ALEFT  AND  ARITE  ARE  END-PTS.  OF  A-INTERVAL  FOR  REGULUS  FALSI  ROUT. 
ALP  IS  THE  TABLE  OF  ATTACK  ANGLES  COMPUTED  IN  THE  OFF-DESIGN  CALC. 
ALPHA  IS  THE  ATTACK  ANCLE  (RADIANS) 

AMDA  IS  A UNIVERSAL  FUNCTION  BOUND  USED  TO  SET  SUBINTERVAL  LENGTHS 
AMU  IS  THE  RATIO  OF  CAVITY  THICKNESSES  AT  XO  AND  THE  TRAILING  EDCE 
AONE  IS  THE  FUNCT.  A1(A)  NEEDED  TO  COMPUTE  SIG1 — SECT5. 

ASQ  IS  A* * 2 . 

ATI  0 IS  A 2 ( A ) REQUIRED  BY  SIC2. 

A 1 IS  THE  CAVITY  PARAMETER  FOR  THE  FIRST  FOIL  (SHOCKLESS  ENTRY). 

A2  IS  THE  CAVITY  PARAMETER  FOR  THE  SECOND  FOIL  (NOSE  SINGULARITY) 
BB  IS  LITTLE  B,  THE  COEFFICIENT  OF  THE  'SQUARE  ROOT'  BEHAVIOR  TERM 
BET  1 AND  BET2  ARE  FUNCTIONS  OF  OMGSQ  AND  A. 

BETA  AS  A VARIABLE  IS  CUNPUTED  IN  S ECT 5 , -MAX ( ( T AU-S IG 2 ) / SI G l ) . 

BINT  IS  THE  PARAMETER  CAPITAL  B. 

BONE  IS  B 1 ( A ) REQUIRED  BY  SIG1. 

BTfO  IS  B 2 ( A ) NEEDED  FOR  SIG2  — SECT5. 

CAY  IS  THE  CAVITATION  NUMBER. 

CD  IS  THE  DRAG  COEFFICIENT. 


CL  IS  THE  LIFT  COEFFICIENT. 


CM  IS  THE 
CO  IS  THE 


DELTA  IS 


A 


MOMENT  COEFFICIENT. 
FIRST  TERM  OF  FOURIER 
SPECIAL  FUNCTION  OF 


EXPANSION. 

THE  CAVITY  LENGTH 


ELL. 


C DELT1  IS  A*BETl/2/ OMGSQ. 

C DELT2  IS  A* BET  2/2 /OMGSQ. 

C DINT  IS  THE  PARAMETER  CAPITAL  D. 

C DLPHA  IS  THE  ATTACK  ANGLE  IN  DEGREES. 

C DROIT  IS  THE  RIGHT-HAND  END-PT.  OF  SUBINTERVAL  COMMUNICATED 

C TO  REGUL  THROUGH  THE  SUBROUTINE  LIST. 

C ECAP  IS  THE  PARAMETER  CAPITAL  E OF  THE  'SlIOCKLESS  ENTRY'  CASE. 

C ELL  IS  THE  CAVITY  LENGTH  = A**2+l. 

C EM  IS  THE  PRESSURE  FUNCTION  MULTIPLIER  COMPUTED  FOR  2ND.  FOIL 

C EPSILN  IS  A SPECIAL  FUNCTION  OF  THE  CAVITY  LENGTH. 

C ETA  1 , ETA 2 , ETA 3 ARE  FUNCTION  VALUES  AT  XC 1 , XC  2 , AND  XC3. 

C ETN 1 IS  THE  CAMBER  FUNCTION  FOR  THE  FIRST  FOIL. 

C ETN  2 IS  THE  CAMBER  FUNCTION  FOR  THE  SECOND  FOIL. 

C ETP1  IS  THE  DERIVATIVE  OF  THE  1ST.  CAMBER  FUNCTION. 

C ETP2  IS  THE  DERIVATIVE  OF  THE  2ND.  CAMBER  FUNCTION. 

C ETPRM  IS  THE  DERIVATIVE  AS  COMPUTED  BY  DIFFER. 

C FX  IS  THE  TABLE  OF  ORDINATES  NEEDED  BY  THE  INTEGRATION  ROUT.  NT CRT 

C G IS  THE  PRESSURE  MINUS  THE  SQ.  ROOT  TERM  - G=P-BB* ( 1 -X ) ** . 5 . 

C GAM  1 , GAM 2 ARE  FUNCTIONS  OF  SMGSQ,  A,  BET  1 , BET  2 . 

C GAUCHE  IS  THE  LEFT-HAND  END  PT . OF  THE  SUBINTERVAL  CONTAINING. 

C ROOT  - SUBINTERVAL  IS  COMMUNICATED  TO  REGUL. 

C H IS  THE  MAPPING  OF  G UNDER  THE  MAPPING  OF  THE  Z-PL.  TO  THE  NU-PL. 

C IFOLSJ  IS  A FLAC  USED  IN  SECT5  - =1,2  FOR  1ST,  2ND  FOILS,  RESP. 

C II  IS  SUBSCRIPT  OF  X DURINC  EVALUATION  OF  FOIL  FUNCTIONS 

C NA  IS  THE  NUMBER  OF  LITTLE  A VALUES-CAVIT Y PARAMETER  VALUF.S-INPUT 

C FOR  THE  OFF  DESIGN  CALCULATION. 

C NM  IS  THE  NO.  OF  ORDINATES  4 HICH  IS  REQUIRED  BY  NTGRTE. 

C NN  IS  THE  NUMBER  OF  MESH  POINTS  IN  INTERVAL  (0,1). 

C NSVJCH=1  FOR  NEXT  CASE  P=CONST.,  =2  FOR  NEXT  CASE  NONCONST.  P, 

C =4  FOR  NO  MORE  DATA 

C OMGSQ  IS  A* SELL 

C ONE  IS  THE  INTEGRAL  OF  C(X)  OVER  (0,1). 

C P IS  THE  INPUT  PRESSURE  FUNCTION  DEFINED  OVER  (0,1),  (Z-PLANE). 

C PCON  IS  THE  VALUE  OF  THE  INPUT  PRESSURE  FOR  A CONST.  P CASE. 

C Q IS  THE  TRANSFORMED  P UNDER  MAPPING  FROM  Z TO  NU  PLANES. 

C SALL  IS  SQ.  ROOT  OF  A * E L L 

C SELL  IS  S>  ROOT  OF  ELL 

C SIC1  IS  THE  FUNCTION  SIGMA1(X,A) — OFF-DESIGN  CALCULATION. 

C SIG2  IS  THE  FUNCTION  SICMA2(X,A) — FOR  THE  OFF-DESICN  CALCULATION. 

C TAU  IS  REDUCED  CAVITY  USED  IN  OFF-DESICN  CALCULATION. 

C TEE  IS  THE  STRENGTH  PARAMETER  APPLIED  AT  THE  TRAILING  EDGE. 

C THETA  IS  THE  CAVITY  CONTOUR  MULTIPLIER  USED  IN  THE  OFF-DESIGN  CALC 

C THR  IS  THE  INTEGRAL  OF  H (T ) OVER(O.l) — IN  NU-PLANE. 

C TI  1 , TI2.TI 3,TI4  ARE  SPECIAL  INTEGRALS  COMPUTED  BY  T I l FN , T I 2 F N , T 1 3F 

C N , TI4FN  FOR  SECTS  FOR  THE  FINITE  CAVITY  CASE. 

C TI5.T16  ARE  SPECIAL  INTEGRALS  COMPUTED  BY  TI5FN.TI6FN  FOR  THE 

C INFINITE  CAVITY  CASE  (CAY=0)  OF  SECT5. 

C T)0  IS  THE  INTEGRAL  OF  X*G(X)  OVER  (0,1). 

C W A IS  THE  FUNCT.  /(A)  COMPUTED  IN  SECT5  — NEEDED  FOR  ALPHA. 

C X IS  THE  INDEPENDENT  VARIABLE.  P IS  INPUT  AS  FUNCT.  OF  X IN  (0,1). 

C XCO,XCl,XC2,XC3  SPECIFY  PTS.  IN  (0,1)  FOR  QUADRATIC  FIT  BY  PARAB 

C XH  IS  THE  STEP-SIZE  — XH*(NN-1)  = 1. 

C XKZ  IS  IMACE  OF  X UNDER  MAPPING  OF  Z-PL.  TO  NU-PL.  — 

C - A* (X/ (ELL-X) )**. 5. 

C XO  IS  THE  THICKNESS  CONTOL  POINT  NEAR  THE  LEADING  EDCE 

C a SO  IS  X * * 2 « 


o o o n 


YYU 1 IS  THE  CAVITY  FUNCTION  FOR  THE  FIRST  FOIL. 

YYU  2 IS  THE  CAVITY  FUNCTION  FOR  THE  SECOND  FOIL. 

INPUT  DATA 

FIRST  CARD  - 

NSNCH  - SET  TO  ONE  (1)  FOR  CONSTANT  PRESSURE  CASE 

- SET  TO  TJO  (2)  FOR  NORMAL  PRESSURE  CASE 
IPUNCH  - ZERO  (0)  - NO  PUNCH 

- ONE  (1)  FOR  PUNCHED  OUTPUT 
K?J  IT  - ZERO  (0)  - RUN  OFF  DESIGN  CONDITIONS 

- ONE  (1)  NO  OFF-DESICN  CALCULATIONS 

IUPPER  - ONE  - UPPER  SURFACE  CONTOUR  IS  TO  BE  CHECKED  FOR 

INTERFERENCE  H ITH  CAVITY 

IPRNT  - ZERO  (0)  PRINTS  OVERALL  DESIGN  PARAMETERS  ONLY. 

- ONE  (I)  PRINTS  DETAILED  RESULTS  LIKE  CAVITY  SHAPE  ETC. 
NEST  - ONE  (1)  FOR  ESTIMATION  MODE 

- ZERO  (0)  FOR  DESIGN  MODE. 

SECOND  CARD 

NN  - NUMBER  OF  INPUT  POINTS  ON  PRESSURE  CURVE 
THIRD  CARD 

XH  - DELTA  X INCREMENTS  ALONG  THE  CHORD  LINE 
AMDA  - LAMDA  FUNCTION  USED  TO  LIMIT  INTEGRATION 
BB  - COEFFICIENT  OF  SQUARE  ROOT  BEHAVIOR 
STEP  - STEP  SIZE  IN  THE  SEARCH  FOR  A 
AMX  - MAXIMUM  VALUE  FOR  A 
FOURTH  CARD 

NA  - NUMBER  OF  OFF-DESIGN  CAVITY  LENGTH  VALUES  TO  BE  COMPUTED 
FIFTH  CARD 

AA ( I ) - ARRAY  FOR  CAVITY  LENGTH  VALUES 
SIXTH  CARD 

ITHET  - NUMBER  OF  THETA  VALUES  TO  BE  COMPUTED 
SEVENTH  CARD 

THET(I)  - ARRAY  FOR  THETA  VALUES 
NEXT  NN/8  CARDS 

SUPPER  (I ) - ARRAY  FOR  COORDINATES  OF  UPPER  SURFACE 
NEXT  NN  CARDS 

X(I)  - ARRAY  CONTAINING  INCREMENTS  OF  CHORD  SPACING 
PP(I)  - ARRAY  OF  PRESSURE  DIAGRAM  POINTS  COORSPONDING  TO  EACH  X 
NEXT  T.J0  CARDS  FOR  AS  MANY  CASES  AS  ARE  TO  BE  RUN 
KAY  - DESIGN  CAVITATION  NUMBER 
CL  - DESIGN  LIFT  COEFFICIENT 

TEE  - MAXIMUN  CAVITY  THICKNESS  AT  TRAILING  EDGE 
S - PEAK  PRESSURE  LOCATION 
NEXT  CARD 

NMU  - NUMBER  OF  MU  VALUES 

AMU V - ARRAY  CONTAINING  THE  MU  VALUES 


REAL  KAY (50 ) 

REAL  LAMDA 
DIMENSION  A M U V ( 5 ) 

DIMENSION  PP(101) 

COMMON  P(101  ) ,C(10l  ) ,Q(101 ) ,H(101 ) ,X(10l ) ,XSQ(101 ) ,FX(10l  ) ,ETN1  (10 
I 1 ) ,ETN2(10I ) ,ETP1 (101 ) ,ETP2(101) ,YYU l (10 l ) ,YYU2( 101 ) ,TAU( 101  ) ,ALP( 
2101 ) .DEE ( 101 ) , AA(50) , KAY, CAY, A, A1 , A 2 , AC  A P, ALPHA , AMDA , ASQ, BR , BINT, C 
30, CL, CM, DELTA, DINT , ECAP , ELL , EM , EPS ILN , ETA  1 , ETA 2 , ET A 3 , TEE , ON E , TI  0 , T 
4HR, THETA, TSAVE,XC0,XC1,XC2,XC3,XH,XKZ,NM,NN, JL6, JM6, JR  6, JL 1 3 , JM13, 
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5JR13, JTTl, JTM1 , JTT2, JTM2, JTT3, JTM3, JTT4, JTM4, JTT5, JTM5, JTT7, JTM7, J 
6TT8, JTM8, JTT9, JTM9, JTTIO, JTMIO, JTT II , JIM  1 1 , JTT 1 2 , JTM 1 2 , JTT 2 1 , JTM21 
7, JTT23, JTM23, JTT24, JTM24, JTT25, JTM25 
COMMON  /BLKl/VIBAR, V1BAR1 ,ECAP1, YC12 
COMMON  /BLK2/LAMDA, S, IPUNCH, IERROR , IPRNT 
COMMON  /BLK3/SEL1 

COMMON  /BLK4/NA, ALPA(8 ) ,THET (8) , ITHET , IUPPER , SUPPER ( 1 0 1 ) 

COMMON  /BLK6/AMU, XO , K , ACL , T , B I,B2,B3,B4,B5,B6,B7,FT,HT,EPSLN, DNM, P 
1 1,  TI  ASQ,SQA,Fll,FlXO,BETA,FE,tIETA,NEND,IXO,NCLM,NEST,NEM 
COMMON  /BLK.8/STEP, AMX 

COMMON  / BLK9 /AMO , AM  1 , IMU  , EMO  , EM  1 , MSN 
C 

IOO  FORMAT  (1012) 

110  FORMAT  (1013) 

120  FORMAT  (5F10. 5) 

130  FORMAT  (2F15. 5) 

140  FORMAT  (36H1  LISTING  OF  INPUT  DATA  FOLLOW S/23H0  NO.  OF  CRID- 

1PTS.  NN-I3.17H  STEP-SIZE  XH=E16.8,24H  FUNCTION  BOUND  A M D A = E 1 

26.8/23H  COEF.  OF  SQ.  ROOT  COMP . =E 1 6 . 3 , 2 7H  CAVITATION  NO.  K ( 

3CAY)=E  16.  8,  16H  LIFT  COEF.  CL=E  1 6 . 8/ 2 1 H PARAMETER  T (TEE)=E16.8, 

45X,'S=',F7.3,5X,'XO=',F7.3) 

150  FORMAT  (8F10.5) 

C SECTION  I.  INPUT  SECTION 

READ  (5,100)  NS'.ICH,  IPUNCH,  KNIT,  IUPPER,  IPRNT,  NEST 

PI-3.  141593 

PCON-l.O 

IF  (NSMCH-3)  1 70,  160,  160 
160  CALL  EXIT 
170  READ  (5,110)  NN.NMU 

READ  (5,120)  XH,AMDA,BB, STEP, AMX 
IF  (KNIT)  180,230,  180 
180  READ  (5,210)  NA , ( AA ( I ) , I- 1 , NA ) 

DO  190  1=1,8 

190  ALP A (I )=FLOAT (I ) 

IF  (IUPPER)  220,200,220 
200  READ  (5,210)  ITHET , (THET (I ), 1=1 , ITHET ) 

CO  TO  230 

210  FORMAT  ( I 5 , / ( 8F  1 0 . 5 ) ) 

220  READ  (5,150)  ( S UP P ER ( I ) , I = 1 , NN  ) 

230  IF  (NS-JCH-1  ) 240,240,260 
240  READ  (5, 120)  PCON 
X(1 )=0. 0 
NNN  =NN - 1 
DO  250  1=1, NNN 

X(I  + 1 )=X(I  )+Xll 
P(I )=PCON 
250  CONTINUE 

P (NN )=PCON 
GO  TO  270 
260  CONTINUE 

270  READ  (5,130)  ( X (I ) , P ( I ) , I = 1 , NN  ) 

READ  (5,120)  KAY(1),CL,TEE,S 
CAY»KAY(1  ) 

LAMDA-0. 05 
NNN-NN-1 
X0=0.  1 
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EPS-1. OE-6 
IMU-1 

FORMAT  (I10.7F10.5) 

MS1!  -1 

IF  (NEST. EQ. 1)  GO  TO  300 
MU-0 

READ  280,  NMU, (AMUV(I) ,I-1,NMU) 

MU-UU+1 

IF  (MU. GT. NMU ) GO  TO  370 
IF  (MSW.EQ.O)  AMU-AMUV (MU ) 

IF  (MS1.EQ.1)  EM-AMUV  ( MU  ) 

NEND-NN 
IXO-1 1 

SECTION  II.  SQUARE  ROOT  BEHAVIOR. 

IF  (NSJCH-1)  310,310,330 
DO  320  1=1, NN 

G (I )-P (I ) 

CONTINUE 
BB-0 . 0 
GO  TO  350 
DO  340  1=1, NN 

H(I)=P(I)/CL 

G(I ) = P (I )-BB*SQRT ( 1 . 0-X(I  ) ) 

WRITE  (6,140)  NN , XH , AMDA , BB , KAY ( 1 ) , CL , TEE , S , XO 
WRITE  (6,150)  ( I , X ( I ) ,P(I)  ,C(I)  ,1  = 1 ,NN) 

IF  (NEST.EQ.l)  GO  TO  360 
NCLM-0 
N EM-0 

CALL  SECT6 
GO  TO  290 
CALL  EST 
CONTINUE 
GO  TO  270 

CALL  SECT4  ( N S 1 CH , PCON , Id  IT ) 

IF  (K!  IT.  EQ.  0.  OR.  IERROR.  EQ.  1 ) GO  TO  270 
CALL  SECTS  (NSJCH,PCON) 

CO  TO  270 
END 


SUBROUTINE  EST 
REAL  KAY (50) , K 

COMMON  P(101) , G ( 1 0 1 ) ,Q(10l) ,H(101) , X ( 1 0 1 ) ,XSQ(101 ) , FX (101) , ETN 1(10 

I 1 ) ,ETN2  (101  ) , ETP1  (101 ) , ETP2 (101 ) , YYU1 (101 ) , YYU2 (1 01 ) , TAU ( 101 ) , ALP( 
2101 ) ,DEE ( 101 ) ,AA( 50) , KAY, CAY, DUM.Al ,A2 , ACAP, ALPHA , AMDA , ASQ, BB , BINT 
3,  CO,  CL  , CM,  DELTA,  DINT  , EC  A P,  ELL  , EM  , E PS  ILN  , ET  A 1 , ET  A 2 , ET  A 3 , TF.  L , ON  E , T!  0 

4 , THR  , THETA  , TS A VE  , XC  0 , XC  1 , XC  2 , XC  3 , XH  , XKZ  , NM  , NN  , JL  6 , JM  6 , JR  6 , JL  1 3 , JM  1 

5 3,  JR13, JTT1 , JTM1 , JTT 2 , JTM2 , JTT 3 , JTM3 , JTT 4 . JTM4 , JTT 5 , JTM5 , JTT 7 , JTM 7 
6,  JTT8,  JTM8,  JTT9,  JTM9,  JTT10,  JT’tlO,  JTT  1 1 , JTM  11  , JTT  1 2 , JTM  1 2 , JTT  2 1 , JTM 
721, JTT 23, JTM23, JTT 24, JTM24, JTT25, JTM25 

COMMON  /BLK6/AMU , XO , K , ACL , T , B 1 , B 2 , B 3 , B 4 , B 5 , B6 , B 7 , FT , HT , EPSLN , DNM , P 

II  ,TJ  ASQ,  SQA,  FI  1 , F 1X0  , BETA,  FE , It  ETA  , NEND  , IXO  , NCLM,  NEST  , NEM 
COMMON  /BLK2/LAMDA , S , I PUNCH, IERROR , 1PRNT 
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COMMON  /BLK9/AM0, AMI ,1MU, EM0.EM1 ,MSN 

IMU-1 

IMF-0 

NEM-0 

NEST-1 

NCLM-1 

CLD-CL 

EM-0 

CALL  SECT6 
CLO-CL 

PRINT  170,  EM, CL 
EM-1 

CALL  SECT6 
CL  1 -CL 

PRINT  170,  EM, CL 

IF  (CL0.LT.CL1)  GO  TO  100 

IMF-1 

TEMP-CLO 

CLO-CL  1 

CL  1 -TEMP 

100  IF  (CLD.GT.CL1)  GO  TO  140 
CL-CLD 
NCLM-0 

PRINT  180,  CLD 

IF  (CLD.GT.CLO)  GO  TO  110 

EM-0 

CALL  SECT6 

PRINT  190,  EM, AMU 

AMO- AMU 

EMO-EM 

EM-1 

CALL  SECT6 
PRINT  190,  EM, AMU 
AMI -AMU 
EMI -EM 
GO  TO  150 

110  IF  (IMF.EQ.l)  GO  TO  120 
EM-1 

CALL  SECT6 
PRINT  190,  EM, AMU 
AMI -AMU 
EMI -EM 
GO  TO  130 
120  EM-0 

CALL  SECT6 
PRINT  190,  EM, AMU 
AMI -AMU 
EMI -EM 
130  NEM-1 

CALL  SECT6 

PRINT  190,  EM, AMU 

AMO-AMU 

EMO-EM 


CO  TO 
140  PRINT 


150 

160 


IMU-0 


n n o 


/ 
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CL-CLD 
EMO-O. 0 
EM  1*0.0 
AMO-O. 0 
AMI-0. 0 
150  CONTINUE 

IF ( IPUNC  H. 
RETURN 

160  FORMAT  (5X 
170  FORMAT  (5X 
180  FORMAT  (5X 
190  FORMAT  (5X 
END 


EQ.l)  PUNCH  1 50 , K , CL , T , S , AMO , AM  1 

DESIGN  CL  IS  GREATER  THAN  MAXIMIUM  PERMISSIBLE  CL 
FOR  M-',E15.7,'  MAXIMUM  VALUE  OF  CL  IS=',E15.7) 

, ' THE  RANGE  OF  M AND  MU  FOR  DESIGN  CL**',E15.7) 

, 'M»' ,E15. 7 , 5X , ' MU=' ,E15. 7) 


) 


SUBROUTINE  SECT6 

MAIN  PROGRAM  FOR  THE  THIRD  DESIGN  PROCEDURE 


REAL  K 
EXTERNAL  FA 

DOUBLE  PRECISION  D A , D ELL , D XO , DQL X , DR F AC , D AMLX , DD E LT A , DE PS LN 
DIMENSION  XOV (3) ,IXOV(3) ,CLV(2) ,AMUV(3) 

REAL  KAY (50) 

COMMON  P ( 1 01 ) , C ( 10 1 ) , Q ( 1 0 1 ) ,H ( 10 1 ) , X( 1 0 1 ) , XSQ ( 10 1 ) , FX ( 1 0 1 ) , ETN  1 ( 10 
11 ) ,ETN2 (101 ) ,ETP1 (101 ) ,ETP2(101 ) ,YYUl (101 ) ,YYH2 (101 ) ,TAU(10l ) ,ALP( 
2101 ) ,DEE ( 101 ) , AA( 50) , KAY, CAY ,DUM, A1 , A 2 , AC A P , ALPH A , AMD A , AS Q , BB , BINT 

3,  CO, CL , CM, DELTA, DINT , EC AP , ELL , EM  , EPS  IL  N , ET  A 1 , ET  A 2 , ET  A 3 , TE  E , ON  E , Tl  0 

4 ,  THR , THETA , TSA VE , XCO , XC 1 , XC  2 , XC  3 , XH , XKZ , NM , NN , JL 6 , JU6 , JR 6 , JL 1 3 , JM 1 
53, JR  13 , JTT 1 , JTM 1 , JTT  2 , JTM2 , JTT  3 , JTM 3 , JTT 4 , JTM4 , JTT  5 , JTM5 , JTT  7 , JTM  7 
6 , JTT 8, JTM 8, JTT 9 , JTM 9 , JTT 1 0 , JTM 1 0 , JTT 1 1 , JTM 1 1 , JTT  1 2 , JTM  1 2 , JTT 2 1 , JTM 
72 1 , JTT 23 , JTM 2 3 , JTT 24 , JTM 24 , JTT  2 5 , JTM 2 5 

COMMON  / BLK6 / AMU , XO , K , ACL , T , B 1 , B 2 , B 3 , B 4 , B 5 , B 6 , B 7 , FT , HT , EPSLN , DNM , P 
II , TJ  ASQ, SQA, Fll , F1X0 , BETA, FE, HETA, NEND , IXO , NCLM , NEST , NEM 
COMMON  /BLK7/DDELTA .DEPSLN 
COMMON  /BLK8/STEP, AMX 

COMMON  /BLK2/LAMDA, S , IPUNCH, IERROR , IPRNT 
COMMON  /3LK1/V1RAR, VI  BARI ,ECAPl,YC12 
COMMON  /BLK3/SEL1 

COMMON  / B L K 9 / AM 0 , AM  1 , IMU  , EMO , EM  1 , MS>f 
DATA  CLV/0.  1,2.0 / 

DATA  XOV/O. 04,0.  1,0.  16/ 

DATA  IXOV/3,6,9/ 

G0(XX)“XX*B3/2  + (SQRT(XX*(ELL-XX))+ELL*ATAN(SQRT (XX/ ( ELL -XX ) ) ) ) *B  1 / 
12. 0-ATAN (SQRT (XX/ ( ELL -XX ) ) ) * B 4+ALO G ( ELL ) * B 6 / 2 . 0- ( BTX 7 (XX ) -BTG 5 ( 1 . 0 
2)/2.0) 

C SECTION  1 . INPUT 

ACL-CL 
T-TEE 

PI-3.  14  1 59265 
K-CAY 

IF  (K. Eg. 0.0)  GO  TO  280 


L > 
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C 2.  ITERATE  TO  GET  THE  VALUE  OF  A 

Al-0.  1 
AL-A1 
AR-Al+1 . 0 
SH-STEP  s 
EPS I-O. 01 

FL-FA (AL  , PCON  , MS'/  CH) 

100  FR-FA(AR, PCON ,NSJ CH) 

TT-FR*FL 

IF  (TT.LT.0.0)  CO  TO  120 

AL-AR 

FL-FR 

AR=AR+SH 

IF  (AR.LT.AMX)  GO  TO  100 
PRINT  110 

110  FORMAT  ( 5X , ' NO  SOLUTION  FOR  A IS  FOUND  IN  THE  GIVEN  RANGE') 

RETURN 

120  CALL  REGUL  ( FA , 0 . 0, AL , AR , 0.  00005, A, PCON, NSW  CH) 

PRINT  130,  A 

130  FORMAT  ('  A=',E15.8) 

C COMPUTE  M,ETC. 

ERR  1 =FA ( A , PCON  , N SW  CH) 

PRINT  140,  ERR  1 

140  FORMAT  ('  ERROR  IN  FA  » '.E15.8) 

IF  (MSI . EQ. 1 ) CO  TO  150 
IF  (NEST. EQ. 1 ) RETURN 

150  BINT=(  1 . 0-Ell)*  CL /PI  / ( 1 . 0+K ) /DNM* ( EPSLN+A*DELTA ) / A/ EL  L / 2 . 0+EM  / A/ ELL 
1/4. 0/PI/C1. 0+K)*B5 

B 2PD= ( 1 . 0-EM)*CL*DELTA/PI/( 1 . O+K) /DNM+EM*B3/4 . 0/PI/( 1 . 0+K) 
DINT=B2P!)-BINT/2. 0 

ECAP^BETA-  ( 1 . 0-EM)  *CL*T.l  ASQ/2. 0/PI  / ( 1 . 0+K ) -EM* B 2 / 2 . 0/PI  /(  1 . 0+K)  / EL 
1L 

C-l. 0/A 
II-NN 

ALPHA=T-( 1-EK) *CL*4 . 0 *A* * 2 *ELL *S QA* F 1 1 /P I / ( 1 . 0+K) /DNM- (A-ELL*ATAK ( 
11.0/A) ) *A*BINT+B2PD-ECAP*A-EM*G0 ( 1.0)/2.0/PI/(l. 0+K) 

ACAP  — 2. 0*( 1 . 0-EM) *CL*SQA/?I / ( 1 . 0+K) /DNM 
CD=2. 0*P1 *( 1 . 0+K ) *ELL*(A*BINT+ECAP/2. 0)**2 
DA=A 

CM— (1. 0-EM  ) * CL  * ( (3. 0D0*DDELTA-4.  OD  0 *D  E PS  LN  *DA  ) +D  A*  * 2 * (DD  ELT  A- 2 . OD 
10*DA*DEPSLN)  ) /8. 0 /D  E PS  LN  / S QP.T  ( ELL  ) -EM*  B I CM  ( A ) 

XBAR  — CM/CL 
YC 1 =A*  ECA  P 
ALOD-CL/CD 

V 1 BAR=ALPUA-B  2PD  + EM*B3/4. 0/PI / ( 1 . 0+K) 

EC AP  1 =EC A P 
SEL1=SQRT (ELL ) 

PRINT  180,  K , AMU , XO 
PRINT  190,  CL  , T , S 

PRINT  200,  ALOD, CM, CD, XBAR, ALPHA, YC 1 , ACAP, EM, BINT, DINT 
PRINT  160,  V1BAR  , F.CAP  1 , SEL  l 

160  FORMAT  ('  INPUT  FOR  OFF-DESIGH  CALCULAT  IONS : ' , / , 5X  , ' V 1 BAR « 

1 ' ,E15.  7. 5X, ' ECAP1-' , E15.  7. 5X, 'SQRT  OF  ELL  1 - ' , E 1 5 . 7 ) 

IF  ( I PUNCH . ECJ.  1 ) PUNCH  170,  K,XO,AMU 

IF  (IPUNCII.  EQ.  1 ) PUNCH  170,  CL  , T , S , ALOD  , CM  , CD  , XBAR  , ALPHA  , ELL  , YC  1 , E 
1M, ACAP, BINT, DINT 
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I FORMAT  (5E15.7) 

I FORMAT  (5X,'K-',F7.4,5X,'AMU=',E15.  7,5X,'1C0-',F7.4) 

I FORMAT  (5X,'CL-' ,E15. 7,5X,'T-',F7.4,5X,'S-' ,F7.4) 

> FORMAT  (5X, 'L/D='  ,E15.7, 5X, 'CM-'  , E 1 5 . 7 , 5X , ' CD - ' , E 1 5 . 7, 5X, 'XBAR-' ,E 
115. 7, 5X, 'ALPHA-' , E 1 5 . 7 , / , 5X , ' YC 1 - ' , E 1 5 . 7 , 5X , ' ACAP- ' ,E15. 7, 5X, 'EM-' 
2,E15.7,5X,'B-',E15.7,5X,'D-',E15.7) 

CONTINUE 

IF  (IPRNT.EQ.O)  GO  TO  270 
COMPUTE  Y,  ETA,  ETC. 

BT5-BTG5 (1.0) 

ETP  1 (1  )=9999. 0 

) FORMAT  (14X, 'X' , I4X, 'Y' , 13X, 'ETA' , 12X, 'ETP' , 12X,  'PTOT'  ) 

PRINT  210 

DO  260  1=1, NN 

FI1-I-1 

FN  1 -NN- 1 

XX-FI1/FN1 

II«(XX-X(1 ) ) /XH+1+0. 5 
SQXX-SQRT (XX*(ELL-XX) ) 

C-SQRT (XX/ (ELL-XX) ) 

XNU=A*C 

TH-ARCOS (1 .0-2. 0*XNU) 

Pl-H(I) 

IF  (I.NE.  1. AND. I.NE.NN)  GO  TO  220 
IF  (I . EQ.  1 ) PM 1-9998. 0 
IF  (I . EQ. NN ) PM  1 =0 . 0 
GO  TO  230 

) PM  1 - 1 . 0 / TAN  (Til  / 2 . 0 ) 

) PTOT=EM*CL*P 1-2 . 0*ACAP*( 1 . 0+K)*PMl 

YYU 1 (I )=ALPHA*XX-2. 0 * A* * 2 *ELL * ACA P * F 1 FUN C ( XN U , A ) +A* B INT * ( S QXX-ELL * 
1ATAN (C) )-3 2PD  *XX+EC AP  *SQXX+EM*G  0 (XX ) / 2 . 0/P1 / ( 1 . 0+K) 

ETN 1 (I )=ALPHA*XX-B2PD *XX+A*  B IWT  * ( ELL  * ATAN ( C ) -S  QXX ) -ECA  P*5  QXX+EM/ 2 . 
10/PI/(  1 . 0+K.)*  (XX*B3/2. 0-(SQXX+ELL*ATAN(C)  ) * B 1 / 2 . 0+ATAN  (C  ) *B4+ALOC  ( 
2 ELL /( ELL-XX )) * CL / 2 . 0- ( ALO C (ELL /SQRT (ELL-XX ))*CL  + FTJX( 1 1, XX )-BT5/2. 
30-BTX7 (XX) ) ) 

IF  (XX. EQ. 0. 0)  GO  TO  240 

ETP1 (I )=ALPHA-BINT*( 1 . 0-2. 0*XNU)/2. 0-D I NT -ECA P * ( E LL -2 .0*XX)/2.0/SQ 
1XX+EM/2 . 0/PI / ( 1 . 0+U) * (B3/2 . 0-B 1 / 2 . 0/C+B4/2 . 0 / SQXX+CL *XX / 2 . 0/(ELL-X 
2X ) +C PVI (I  I , XX) /2 . 0 /SQRT (XX ) / ( ELL -XX ) ) 

) IF  (IPRNT.EQ.l)  PRINT  250  , 1 1 , XX , Y Y U 1 ( I ) , ET N 1 ( I ) , ET P 1 ( I ) , PTOT 

IF  (IPUNCH.  EQ.  1)  PUNCH  170,  XX , YYU  1 ( I ) , ETN 1 ( I ) , ET P 1 ( I ) , PTOT 
) FORMAT  ( 5X , 13, 5E15.  7) 

) CONTINUE 
) CONTINUE 
RETURN 

CASE  K-0 
) CONTINUE 

DO  290  1=1, NN 

P (I ) -H (I ) 

) C(I)-P(I) 

STAR-0.66 
BT5-BTC5 (1.0) 

BT8-BTG8 (1.0) 

BT9-BTG9 ( 1 . 0) 

FT11-FTJX(NEND,  1. 0) 

B11-FT11-BT5/2.0-BT8 
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B12«2.0*BT8-FTll 
SQ2-SQRT (2. 0) 

IF  (NCLM.EQ.O)  GO  TO  300 

CLTM-(1. 0— EM ) * ( 3. 0*SQ2-AL0G ( 1 . 0+SQ2 ) ) /PI+EM/2. 0/PI*(2.0*BT9-B12) 

CLTM-1. 0/CLTM 

CL=CLTM*T 

IF  (EM.EQ. I.  . AND. S. GT. STAR)  CL-99.9 
RETURN 

300  IF  (NEM.EQ.O)  CO  TO  310 

TRSQ=3. 0*SQ2-ALOG(1 . 0+SQ2) 

EM-(TRSQ-PI*T/CL) /(TRSQ-BT9+B  12/2. 0) 

310  CLOT-CL/T 

sqxo=sqRT (xo ) 

B1X0=FTJX(IX0 , XO )-BT5/2 . 0-BTG8 (SqXO ) 

cn-2.o*((3. 0*sq2-ALoc ( l . o+sq2) ) *(o.  5+o.  5/sqxo ) - (i . o+2.  o*sqxo )*SqR 
IT (1 . o+sqxo ) /X0**0. 75+ALOC (X0**0. 25+SqRT ( 1 . O+SqXO ) ) /XO ) 
CI»BT8/SqXO-BTG8(SqXO) /X0-(0. 5-0. 5 /XO ) * BT 5+ ( 0 . 5-0. 5/SqXO)* FT  11 
CJ-C11+CI 

IF  (MSI . Eq. 0)  GO  TO  320 

AMU«X0*(EM*CJ-C11 ) /2. 0 /P I *CLOT+XO / 2 . O+SqXO /2. 0 
GO  TO  340 

320  IF  (NEST.Eq.O)  GO  TO  330 

AMU*XO * ( EM*  C J-C 1 1 ) / 2 . 0/PI *CLOT+XO /2 . O+SqXO /2. 0 
RETURN 

330  EM= ( C 1 1 +T * 2 . 0 *P I /CL * ( AMU /XO -0 . 5-0. 5/SqXO) ) /CJ 

340  YC1T-1 . 0-(  1 . 0-EII)  *CLOT/PI*(  3. 0*Sq2-ALOG  ( 1 . 04  S q 2 ) ) -EM*  BT9  *CLOT  / PI+E 
lM*B12*CLOT/2. 0/PI 
YClT*=YClT/2. 0 
YC  1 “YC  IT  *T 
BT5=BT5  *CL 
BT8*=BT  3*CL 
BT9=BT9*CL 
FT11=FT11 *CL 
B11«=B11*CL 
B 1 2 =B  12  *CL 
DO  350  1=1, NN 

P(I)=H(I)*CL 
350  G(I)  = P(I  ) 

ALPHA  = 2 . 0*(1. 0-EM  )*CL/PI+YCl+EM*BT9/PI/2. O+EM/2. 0 /P I * ( FT  1 1 -B T8-BT 5 
1/2.0) 

CD  = ( 2 . 0 * ( 1 . 0-EM)*CL  + EM*BT9/2. 0+P I *YC  1 ) * * 2 /P  I / 2 . 0 
ACAP— 2. 0*(  1.  0-EM )*CL/PI/(  1 . 0+K) 

ALOD  =C  L /CD 

CM— (1. 0-EM)* CL *5. 0/16. 0-EM*BICM(  1.0) 

XBAR  — CM/CL 

V1BAR-ALPHA-2 , 0*(1. 0-EM)*CL/PI 
EC A P 1 =Y  C 1 
S EL  1 = 1 . 0 

PP.INT  180,  K , AMU  , XO 
PRINT  190,  CL , T , S 

PRINT  360,  ALOD, CM, CD, XBAR, ALPHA, YC1, ACAP, EM 
PRINT  160,  V1RAR, ECAP1  ,SELi 

3 60  FORMAT  ( 5X , 'L /D  , E 1 5 . 7 , 5X , ' CM- ' , E 1 5 . 7 , 5X  , ' CD-  ' , E 1 5 . 7, 5X,  'XBAP.-'  ,E 
115.7,  5X,  'ALPHA  = ' , E 1 5 . 7 , / , 5X  , ' YC  1 = ' , F.  1 5 . 7 , 5X,  'ACAP-'  , E 1 5 . 7 , 5X  , ' EM  = ' 
2 , E 1 5 . 7) 

IF  (IPUNCH. Eq.  1 ) PUNCH  170,  K.XO.AMU 
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IF  (IPUNCH. EQ.  1)  PUNCH  170,  CL, T , S , ALOD , CM, CD , XBAR , ALPHA , YC 1 , EM, AC 
1AP 

IF  (IPRNT.EQ.O)  GO  TO  410 

ETPl  (D-9999. 0 

PRINT  210 

DO  400  I-l,NN 

FI  1*1-1 

FNl-NH-1 

XX-FI1/FN1 

II-(XX-X(1 )) /XH+1+0.  5 
SQX-SQRT (XX) 

BT8X-BTC8 (SQX) 

XNU-SQX 

TH-ARCOS ( 1 . 0-2. 0*XN U) 

Pl-H(I) 

IF  (I.NE.  l.AND. I.NE.NN)  GO  TO  370 
IF  (I. EQ.  1)  PM  1=9998. 0 
IF  (I.EQ.NN)  PM1-0.0 
CO  TO  380 

370  PM  1 = 1.0 /TAN (TH / 2 . 0 ) 

380  PT0T=EM*CL*Pl-2. 0*ACAP*(1  0+K)*PMl 

YYU 1 (I  )«ALPHA*XX+( 1 . 0-EM ) * CL /P I * ( ( 1 . 0+2. 0 *S QX ) *XX ** 0 . 2 5*SQRT ( 1 . 0+S 
lQX)-2. 0*XX— ALOG (XX**0. 2 5+SQRT (1 . 0+SQX) ) ) +YC 1 * S QX+EM/ 2 . 0/PI*(SQX*BT 
29+BT5/2 . 0-BT8X) 

ETN  1 (I  )=  (ALPIIA-2 . 0*(  1 . 0-EM)  *CL/PI ) *XX-YC  1 *S  QX-EM*  S QX  / 2 . 0/PI*BT9-EM 
1/2. 0 /PI  * ( FT JX ( 1 1 , XX ) - BT5 /2 . 0-BT8X) 

IF  (XX. EQ. 0. 0)  GO  TO  390 

ETP1  (I) -ALPIIA-2. 0*(  1 . 0-EM)  * CL /PI -YC  1/2. 0/SQX+EM/4. 0/PI  /SQX*  ( CP  VI  (I 
II ,XX)-BT9) 

390  IF  (IPRNT.EQ.l)  PRINT  250  , II  , XX , YYU 1 (I ) , ETN 1 ( I ) , ET P 1 ( I ) , PTOT 
IF  (IPUNCH. EQ.  1)  PUNCH  170,  XX , YYU 1 ( I ) , ETN  1 ( I ) , ET P L ( I ) , PTOT 
400  CONTINUE 
410  RETURN 
END 


FUNCTION  FA  ( A , PCON  , NS'l  C H ) 

DOUBLE  PRECISION  DELL , DTRM 1 , DTRM 2 , DD ELT A , DE PSLN , DK , DB ET A , DA , DXO , DQ 
1LX,DBFAC,DAMLX 
REAL  K 
REAL  KAY (50) 

COMMON  P(101)  , C ( 1 0 1 ) , Q ( 1 0 1 ) ,11  (101)  , X ( 1 0 1 ) ,XSQ(101)  , FX  ( 1 0 1 ) , ETN  1(10 
1 1 ) ,ETN2 ( 101 ) , ETP1 (101 ) ,ETP2(10l ) , YYU 1 (101 ) , YYU 2 (101 ) ,TAU( 101 ) ,ALP( 
2101 ) ,DEE ( 101 ) , AA(50) , KAY , CAY , DUM , A 1 , A2 , AC AP, ALPHA , AMD A, ASQ, BB , BINT 

3,  CO, CL, CM, DELTA , DINT , EC A P, ELL , EM , EPS  I LN , ETA  1 , ET  A 2 , ET  A 3 , TE  E , ON  E , T J 0 

4 ,  THR , THETA , TSAVE , XCO , XC 1 , XC 2 , XC  3 , XH , XKZ , NM , NN , JL  6 , JM 6 , JR 6 , JL 1 3 , JM 1 
5 3, JR  13 , JTT 1 , JTM1 , JTT  2 , JTM  2 , JTT  3 , JTM  3 , JTT4, JTM4 , JTT5, JTM5 , JTT  7 , JTM7 
6, JTT8, JTM8, JTT9, JTM9, JTTIO, JTK10, JTT11 , JTM 1 1 , JTT 1 2 , JTM 1 2 , JTT 2 l , JTM 
7 2 1 , JTT  2 3 , JTM 2 3 , JTT 24 , JTM 2 4 , JTT 2 5 , JTM 2 5 

COMMON  /BLK6/AMU , XO , K , ACL , T , B 1 , B 2 , B 3 , B 4 , B 5 , B 6 , B 7 , FT , HT , EPSLN , DNM , P 
1 1 , TJ  ASQ  , SQA  , F 1 1 , F 1X0  , BETA  , FE,  NET  A , NF.ND  , IXO  , NCLM  , NEST  , NEM 
COMMON  /BLK7/D DELTA , D EPSLN 
COMMON  /BLK9/AM0, AMI , IUU , EMO, EMI  ,MSJ 
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G0(XX)-XX*B3/2+(SQRT (XX* ( ELL-XX ) ) +ELL * ATAN ( S QRT (XX/ ( ELL-XX ) ) ) ) * B 1 / 
12. 0-ATAH  (SQRT  (XX/ ( ELL-XX  ) ) ) *B  4+ALOG  ( ELL  ) * B 6 / 2 . 0-(BTX7  (XX)-BTG5  ( 1 . 0 
2 ) / 2 . 0 ) 

PI-3. 14159265 
SQA-SQRT (A) 

DK-K 

DBETA-1 . 0-1 . 0/DSQRT (1 . 0+DK) 

BETA-DBETA 
ELL- ( 1 . 0+A*  * 2 ) 

DELL-ELL 

DTRMi-DSQRT (DSQRT (DELLJ+l.  0) 

DTRM2-DSQRT (DSQRT (DELL)- l . 0) 

DDELTA-DTRM 1+DIRM2 

DEPSLN-DTRM1-DTRM2 

DN11«DDELTA-A*DEPSLN 

DELTA-DDELTA 

EPSLN-DEPSLN 

DO  100  1=1, NN 

P (I ) -H ( I ) 

100  G ( I ) -P ( I ) 

IF  (NCLM.  EQ.  0)  GO  110 

CLDNM- ( 1 . 0-EM)* ( (2. 0*A**2  + 3. 0)*EPSLN  +A*D  ELT  A ) /DNH+EM* B TG 1 ( A ) 
CL-BETA*  2 . 0*PI*( 1 . 0+K ) *ELL /CLDNM 
CO  TO  120 

110  IF  (MEM. EQ. 0 ) CO  TO  120 

TSQ=(  (2. 0*A*  * 2 + 3 ) * E PS  LM+A*D  ELT  A ) /DM11 

EM- (BETA* 2. 0*PI*( 1 . 0+K ) * ELL / CL-T S Q ) / ( BTC  1 ( A ) -T S Q ) 

120  DO  130  1=1, NN 

P(I)-H(I)*CL 
130  G ( I ) -P (I ) 

DUM=A 

C COMPUTE  THE  VARIOUS  INTEGRALS 

B 1-BINT  1 (DUM) 

B2-BTG1 (A) 

B3-BINT3 (DUM) 

B4-BINT4 (DUM) 

B5  = ELL*B  1-2. 0*B4 
B6-CL 

B7-ALOC (ELL) *B6+FT JX(NEND ,1.0)-2.0*BTX 7(1.0) 

C COMPUTE  H (ETA) 

F11-F1FUNC ( 1 . 0, A) 

HETA-CL* ( 4 . 0*A**2*ELL*SOA*F1  1 /DNM+  ( A-ELL  *ATAN  ( 1 . 0/A)  )*(EPSLN+A*DEL 
1TA) /ELL/ONM-A* ( (2. 0 *A* * 2+ 3 ) * E PS LN+A *D ELTA ) / EL L /D NM ) - ( A-ELL * ATAN ( 1 / 
2A))*B5/2. 0 / ELL+A*  3 2/  ELL-(A+ELL*ATAN ( 1 . 0 / A ) ) * B 1 / 2 . 0+ATAN ( 1 . 0/A)*B4- 
3B7/2. 0 

FE-T-2 . 0*A*SETA+CL* ( A* ( ( 2. 0* A* * 2+3 ) * E PS LN+A*D ELT A ) /ELL/DNM-4 . 0*A** 
1 2*ELL*SQA*F  1 1 / DM  M-  (A-E1.L*  ATAN  ( 1 . 0/A)  ) * ( E PS  LN  +A*D  ELT  A ) / ELL  /DNM ) /PI  / 
2(1. 0+K) 

HE  1 -CL*4 . 0*A**2*ELL*SQA*F1  1 /DNM 

HE  2- (A-ELL *AT AN ( 1 . 0 / A ) ) * CL  * ( E PS LN +A* D ELT A ) / ELL /DNM 
1IE3—  (A-ELL*ATAN(1 . 0/A)  )*B5/2. 0/ELL 
HE4--A*CL*( (2. 0*A**2  + 3)*EPSLN  + A*  D ELI A ) / ELL /DNM 
HE5-A*P<  2/ELL 

H E 6 »- ( A+ELL * AT  AN (1.0/A))*Bt/2.0 
HE  7-ATAN ( 1 . 0 / A ) * 3 4 
HE8--B7/2. 0 
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H9-FTJX(NEND,  1.0) 

H 10«BTX7 (1.0) 

C PRINT  100,HE1,HE2,HE3,HE4,HE5,HE6,HE7,HE8,H9,H10 

H ETA “HE 1+HE2+HE3+HE4+HE5+HE6+HE  7+HE8 
IF  (NCLM.EQ.O)  CO  TO  140 
FA-FE+EM*HETA/PI/ ( 1+K) 

RETURN 

C COMPUTE  THE  FUNCTIONS  F 1 ( 1 , A) , F 1 (X , A) 

140  XNU«A*SQRT (XO/(ELL-XO) ) 

FIXO-FIFUNC(XNU.A) 

II-NN 
XX-1 . 0 
G1A»C0(1. 0) 

II-IXO 

XX«=XO 

GXOA=GO (XO ) 

SQLX=SQRT ( ( ELL-XO ) /XO ) 

Ti ASQ=( (2. 0*A**2+3)*EPSLN+A*DELTA) /ELL/DNM 
DA=A 

DELL=DA**2+1 . 0 
DXO  =XO 

DQLX=DSQRT ( (DELL-DXO ) /DXO ) 

DBFAC=DA-DQLX-DELL*DATAN  ( 1 . 0 /D  A ) +D  E LL  *T)  ATAN  ( 1 . 0 /D  QL  X ) /D  XO 
DAMLX=DA-DQL  X 

FT=T-CL*4.  0*A**2*ELL*SQA* (FI 1-F  1X0 /XO  ) /PI  / ( 1 . O+IC  ) /DNM-DPjF  AC  * CL  * ( E P 
1 SLN+A*DELT A ) / 2 . 0/PI /( 1 . 0+K ) / ELL /DNM- ( B ETA-CL *T1 ASQ/2.0/PI/(l. 0+K.) ) 
2 *DAMLX 

HT  =CL  *4 . 0*A**2*ELL*SQA* (FI 1-F 1X0 /XO ) /DNM+DBF AC*(CL*( E PS  LN+A*DELT  A ) 
1/2. 0/DNM/ELL-B5/4. 0/ELL)  + (-CL  *T-4  ASQ/2. 0+B2/2 . 0 /ELL  ) *D AMLX-  (C  1A-CXO 
2 A/XO  ) / 2 . 0 

IF  (MSI  . EQ.  1 ) GO  TO  1 50 
IF  (NEST. EQ. 0)  GO  TO  160 
150  FA=FE+EM*H ETA/PI / ( 1+K) 

AMU=XO  * ( FT+EM*HT / PI / ( 1+K ) ) /T 
RETURN 

160  FT=FT-AMU*T/XO 

EM=-FT*P I * ( 1 +K ) / HT 
FA=FE+EM*HETA/PI / ( 1+K) 

RETURN 

END 


1 


I 


FUNCTION  CPVI  (II, XX) 

C 

C 

C 

C THE  CAUCHY  PRINCIPAL  VALUE  INTEGRAL  (NU-PLANF.)  APPEARING  IN  THE 

C EQUATION  FOR  THE  DERIVATIVE  IS  TRANSFORMED  TO  THE  SUM  OF  A J F.LL- 

C BEHAVED  INTEGRAL  (DIFF1)  AND  ANOTHER  CAUCHY  PRINCIPAL  VALUE 

C INTEGRAL  J HEN  THE  NU-PLANE  IS  MAPPED  INTO  THF.  Z-PLANE.  THIS 

C ROUTINE  EVALUATES  THE  LATTER  INTEGRAL. 


o o 


NEEDS  COMMON. 

REAL  KAY (50) 

COMMON  P(101)  ,G(101) , Q ( 1 0 l ) ,H(101) , X ( 1 0 1 ) ,XSQ(101 ) ,FX(101)  , ETN  1(10 
11) , ETN2C101 ) ,ETP1 (101 ) ,ETP2(101 ) , YYU 1 ( 1 0 1 ) , YY U 2 ( 1 0 1 ) ,TAU (101 ) ,ALP( 
2101) ,DEE(101) , A A (50) , KAY , CAY , A , A1 , A2 , ACAP , ALPHA , AMDA , ASQ , B B , BINT , C 
30 , CL , CM , DELT A , DINT , EC AP , ELL , EM , EPS ILN , ET A 1 , ET A2 , ET A3 , TE E , ON E , TW  0 , T 
4HR  , THETA , TS A VE , XC  0 , XC 1 , XC  2 , XC  3 , XH , XKZ  , NM , NN , JL  6 , JM6 , JR  6 , JL  1 3 , JM 1 3 , 
5JR13, JTT1, JTM1, JTT2, JTM2, JTT3, JTM3, JTT4, JTM4, JTT5, JTU5, JTT7, JTM7, J 
6TT8, JTM8, JTT9, JTM9, JTT10, JTM10, JTT 1 1 , JTMl 1 , JTT  12  , JTM12 , JTT21 , JTM21 
7 , JTT 23 , JTM23 , JTT 24 , JTM24 , JTT 2 5 , JTM25 

FUN  (XL  , XR , QO , Q 1 , Q2 ) = Q0  *CKLOG ( (XR-XX ) / (XL-XX ) ) +Q 1 * (XR-XL )+Q2  * ( (XR-X 
lX)**2-(XL-XX)**2)/2. 0 
NN2-NN-2 

TEMPI -0. 4142136/2.  42136 

IF  (II.CT.NN)  GO  T 120 
IF  ( I I- 1 ) 120,100,130 
100  CALL  D IFF  5 (XX , 1 , P , CAY , AAO , AA 1 , AA2 ) 

TD1=DIFF2(XX, 3,NN2,P,CAY) 

TD2-DIFF2 (XX, NN2,NN , C.CAY) 

CALL  DIFF8  (XX , CAY , BBO , BB  1 , BB 2 ) 

CPVI=TD1+TD2-AAO*CKLOG((X(1 ) -XX ) / ( X ( 3 ) -XX ) ) +A A 1 * 2 . 0*XH+AA2*( (X(3)- 
1 XX )**2—(X(l ) -XX ) * *2 ) /2 . O+B  B*DIFF4 (X ( 1 ) , XH , B B 0 , BB 1 , 3B 2 ) 

110  RETURN 

120  CPVI-DIFF2 (XX, 1 ,NN,P,CAY) 

CO  TO  110 

130  IF  ( 1 1 -2  ) 100,  140,  150 

140  CALL  DIFF5  (XX , 1 , P , CAY , AAO , AA 1 , AA 2 ) 

TD 1=0  IFF  2 (XX, 3.NN2, P.CAY) 

T D 2 =D  I FF  2 ( X X , NN  2 , NN  , C , C A Y ) 

CALL  DIFF8  (XX , CAY , BBO , BR  1 , BB 2 ) 

CPVI-TD 1+TD2+FUN (X ( 1 ) , X ( 3 ) , AAO , AA 1 , AA 2 )+B 3*D I FF 4 (X ( 1 ) , XH , BB 0 , BB 1 , B 
132) 

GO  TO  110 

150  IF  (II -NN  +4 ) 160,  160,  1 70 
160  XC0=XX 

II M 1=11-1 
IIP1-II+1 

TD1=DIFF2(XX, 1 , IIM1 , P.CAY) 

CALL  DIFF5  (XX , IIK 1 , P , CAY , BBO , BB  1 , BB 2 ) 

TD2=DIFF2(XX,  II ? 1 , NN 2 , P , CAY ) 

TD3=DIFF2( XX , NN  2 , NN , C.CAY) 

CALL  DIFF8  (XX , CAY , CCO , CC  1 , CC 2 ) 

CPVI=TD1+TD2+TD 3+FUN ( X ( I IM 1 ) , X ( I I P 1 ) , BB 0 , BB 1 , BB 2 ) +BB* D I FF 4 (XC 0 , XH , 
1CC0.CC1.CC2) 

CO  TO  110 

170  IF  (II— (NN  — 3) ) 1 60,  180  , 190 
180  II M 1=11-1 
I1P1-II+1 
XCO-O. 0 

TD 1 =D I F F 2 (XX , 1 , IIM1 , P.CAY) 

CALL  DIFF5  ( XX , I I M 1 , P , C A Y , BB 0 , BB 1 , BB 2 ) 

TD  2=DIFF2 (XX , NN  2 , NN , C , CAY ) 

CALL  DIFF8  ( XX , CAY , CCO , CC  1 ,CC2) 

CPVI-TP1+TD2  + FUN (X ( I IM 1 ) , X ( 1 1 P 1 ) , BB 0 , BR 1 , BB 2 ) +BB *D IFF 4 (XX , XH , CCO , C 


200  XCO-O.O 
I IM 1*11-1 
IIP1-II+1 

TD1=DIFF2(XX, 1 , 1 IM  1 , P , CAY ) 

CALL  D IFF  5 (XX , IIM  1 , P , CAY ,B30,BBl,BB2) 

CALL  D IFF  5 (XX , NN 2 , G , CAY , CC 0 , CC  1 , CC 2 ) 

CALL  D IFF  8 (XX , CAY , DD 0 , DD  1 , DD 2 ) 

CPVI-TD1+FUN (X(IIMl) , X (I  I P I ) , BB 0 , BB 1 , BB 2 )+F UM (X (NN - 1 ) , X (NN ) , CC 0 , CC 
1 1 ,CC2)+(-DD0*SQRT  (2. 0*XII ) * ( ALOG  (TEMPI  ) + 1 . 4 142  136) + 2. 0*DD  1 *XH *XH **0 
2. 5/3. 0+14. 0*DD2*Xll**2.  5/15.0)*3B 
GO  TO  110 

210  IF  (II-(NN-l))  200,220,230 
220  XC0=0. 0 
IIM  1 = 11-1 
IIP1-II+1 

TD 1 =D IF  F 2 (XX , 1 , 1 IM 1 , P , CAY ) 

CALL  DIFF5  (XX , IIM  1 , G , CAY , BB 0 , SB  1 , BB 2 ) 

CALL  DIFF8  (XX , CAY , CCO , CC 1 , CC 2 ) 

CPVI-TD1+FUH (X(IIM1) ,X(IIP1) , BB  0 , BB 1 , BB  2 ) +BB*D IFF  4 (XX , XH , CC  0 , CC 1 , C 
1C  2 ) 

CO  TO  110 
230  XCO-O.O 

TD  1 =D IFF  2 (XX , 1 ,NN2,P,CAY) 

XN “X (NN  ) 

XN2=X(NN-2) 

CALL  DIFF5  (XX , NN 2 , G , CA Y , BB 0 , 3B 1 , BB 2 ) 

CALL  DIFF8  (XX , CAY , CC 0 , CC 1 , CC 2 ) 

CPVI =TD 1 +FUN (XN2, XN,BB0, BB1 , BB 2 ) +3 B *D I FF 4 ( 1 . 0,XH,CC0,CC 1 ,CC2) 

GO  TO  110 
END 


FUNCTION  DIFF2  (XX , ILL , IRR , F , CAY ) 

C D I F 2 

C 
C 
C 

C EVALUATE  A DEF.  INTEGRAL  0 VE R (X ( I LL ) , X ( IRR ) WHOSE  INTEGRAND 

C IS  C(X(II) ,U)*F(U)/(U-X(II))  WHERE  U IS  THE  DUMMY 

C VARIABLE,  F IS  P OR  G,  AND  G IS  GIVEN  BY  DIFF7. — NEEDED  BY  ZZZZ. 

C 

C 

C 

C HEEDS  COMMON 

REAL  KAY ( 50 ) 

COMMON  P ( 1 0 1 ) ,G(10l ) , Q ( 1 0 1 ) ,H (101 ) ,X(10  1 ) ,XSQ( 101 ) ,FX( 101 ) ,ETN1  (10 
1 1 ) , ETN  2(101) , ETP1(101),ETP  2(101) ,YYU 1 ( l 0 l ) , YYU 2 ( 1 0 1 ) , TAU ( 1 0 1 ) , ALP ( 
2 101 ) ,DEE ( 10 l ) , AA( 50)  , KAY , DUM , A , A 1 , A2 , ACAP , ALPHA , AM DA , ASQ , BB , BINT , C 
30  , CL,  CM,  DELTA,  DINT  , ECAP,  ELI, , EM,  EPSILN  , ETAl  , ET  A 2 , ET  A 3 , TEE,  TGI  , TG  2 , T 
4G  3 , THETA , TSA  VE  , XC  0 , XC  1 , XC  2 , XC  3 , XH  , XKZ  , MM , NN  , JL6  , JM6  , JR 6 , JL  1 3 , JM  1 3 , 
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5 JR  1 3 , JTT 1 , JTM I , JTT 2 , JTM2 , JTT 3 , JTM3 , JTT4 , JTM  4 , JTT 5 , JTM 5 s JTT  7 , JTM 7 , J 
6TT8,  JTM8,  JTT9,  JTM9,  JTT  10,  JTH10,  JTT  11  , JTM  11  , JTT  1 2 , JTM  1 2 , JTT  21  , JTM  21 
7, JTT 23, JTM23, JTT 24, JTM24 , JTT25, JTM25 
DIMENSION  F(101) 

NM-IRR-ILL+1 

J-ILL 

DO  100  I-l.NM 

FX(I )-DIFF7(XX,X(J ) , ELL , C AY ) * F ( J ) / ( X ( J ) -XX ) 

J-J  + l 

100  CONTINUE 

CALL  NTGRTE  (DIFF2) 

RETURN 

END 


FUNCTION  DIFF4  (XX , XH , BB 0 , BB 1 , BB 2 ) 
D I F 4 


EVALUATES  THE  EXPRESSION  OBTAINED  FOR  THE  INTEGRAL  OVER  (1-2H.1 
OF  F(XX,U)*SQP.T  (1-U)/(U-XX)/SQRT  (U)  WHEN  F/SQRT(U)  IS  REPLACED  BY 
BBO+BB 1*(U-XCO)+332*(U-XCO)**2  OBTAINED  USING  PARAB  PREVIOUSLY. 

(F  IS  GIVEN  BY  DIFF7). 


Tl-SQRT ( 1 . 0-XX) 

T 2“SQRT  ( 2 . 0*XI1) 

DIFF4=BB0*(-T 1 *CKLOG ( (T 1 -T 2 ) /(T1+T2) )-2. 0*T2)+2.0*BB1*T2*(2. 0*XH) / 
13. 0+2. 0*BB2*( ( 1 . 0-XX ) *T  2 * ( 2 . 0*XH)/3. 0-T2*(2. 0*XH)**2/5. 0) 

RETURN 

END 


SUBROUTINE  DIFF5  (XX , ILL , F , CAY , AAO , AA 1 , AA2 ) 

D IF  5 

USES  PARAB  TO  APPROXIMATE  2 ( XX , U ) * F ( U ) / SQRT (U ) BY 
A 0+A  1 * ( U -XX  ) +A  2 * ( U -XX  ) * * 2 
ON  (U(ILL)  ,U(ILL  + 2))~ F-P  OR  G,  Z(XX,U)  IS  DIFF7. 


HEEDS  COMMON 
REAL  KAY ( 50 ) 

COMMON  P,G,O,H,X,XS0,FX,F.TNl  , ETN  2 , ET  P 1 , ET  P 2 , YYU  1 . YYU  2 . TAU  , ALP  , DEE  , 
1AA,  KAY,  DUM,  A,  A I , A 2 , AC  A P , AL  Pll  A , AMDA  , AS  Q . B B , B I NT  , CO , CL , CM  , DELTA , DI  NT 
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2.ECAP, ELL, EM, EPSILN, ETAl , ETA2 , ETA3 , TEE , TC 1 , TC 2 , TC 3 , TH ET A , TS A VE , XC 0 

3 ,XC  1 ,XC2  , XC3,  Xll,  XKZ  , NM,  NN  , JL  6 , JM6  , JR  6 , JL  1 3 , JM  1 3 , JR  1 3 , JTT  I , JTM1  , JTT 

4 2,  JTII2,  JTT3,  JTM3,  JTT4,  JTM4,  JTT5,  JTM5,  JTT7,  JTM7,  JTT8,  JTM8.  JTT9,  JTM9 

5,  JTT 10, JTM  10, JTT  1 1 , JTM1 1 , JTT 1 2 , JTM 1 2 , JTT 2 l , JTM21 , JTT 2 3 , JTM23 , JTT  24 

6 ,  JTM24 , JTT  2 5 , JTM25 

DIMENSION  P(L01) , G ( 1 0 1 ) , Q ( 1 0 1 ) , H ( 1 0 1 ) , X ( 1 0 1 ) , XSQ(lOl) , FX ( 1 0 1 ) , ETN 1 
1(101) ,ETN2( 101) ,ETP 1(101 ) ,ETP2( 101) ,YYU 1(101) , YYU 2 ( 1 0 1 ) , TAU ( 1 0 1 ) ,A 
2LP (101) ,DEE(10l), A A (50) 

DIMENSION  F(101) 

C 

xco-xx 

XC  1**X  (ILL  ) 

XC  2*X ( ILL  + 1 ) 

XC  3“X ( I L L + 2 ) 

ETA  1 *D IFF  7 (XX , XC 1 , ELL , CAY ) * F ( I LL ) 

ETA2*=DIFF7 (XX, XC2, ELL ,CAY)*F(ILL+l ) 

ETA3“DIFF7(XX,XC3,ELL, CAY )*F(ILL+2) 

CALL  PARAB  ( AAO , AA 1 , AA2 ) 

RETURN 

END 


FUNCTION  DIFF7  (XX , UU , ELL , CAY ) 

C EVAUATES  THE  FUNCTION  Z(XX,UU)= 

C SQRT (UU * ( ELL -XX ) ) +S QRT (XX* ( ELL-UU ) ) , K NONZERO, 

C AND  SQRT (UU ) + SQRT (XX)  , K-0. 

C 

c 

C NO  COMMON 

c 

IF  (CAY)  100,120,100 

100  D IFF  7 “SQRT (UU*(ELL-XX) )+SQRT (XX* ( ELL-UU ) ) 

DIFF  7“D I F F 7 *S QRT (ELL-UU ) 

110  RETURN 

120  DIFF7-SQRT (UU )+SQRT (XX) 

CO  TO  110 
END 


C 

C 

C 

C 

C 

C 

C 


SUBROUTINE  DIFF8  (XX , CAY , AAO , AA 1 , AA 2 ) 
DIFS 


APPROXIMATES  Z(XX,U)  BY  AO+A 1 (U -XX )+A 2 (U-XX ) * * 2 OVER 

( 1-2 H,  1 )--Z (X , U)=DIFF7 (X, U , ELL, CAY)  . USES  PARAR. 


no  000000000000  o 
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NEEDS  COMMON 
REAL  KAY ( 50 ) 

COMMON  P ( 1 0 1 ) , G ( 1 0 1 ) , Q ( 1 0 1 ) , H ( 1 0 1 ) , X ( 1 0 1 ) ,XSQ(101) ,FX(101)  , ETN  1 (10 
II) ,ETN2( 101)  ,ETP  1(101  ) ,ETP2(101) ,YYU 1(101) ,YYU2(101) ,TAU(101) ,ALP( 
2101),  DEE  (101),  AA  (50)  , KAY  , DUM  , A , A 1 , A2  , ACAP  , ALPHA  , AMDA  , ASQ  , BB  , B ItiT  , C 
30, CL, CM, DELTA, DINT , EC AP , ELL , EM , EPS ILN , ET A 1 , ET A2 , ET A 3 , TE E , TC  1 , TG 2 , T 
4G  3 , THETA  , TSA  VE  , XC 0 , XC  1 , XC  2 , XC  3 , XII , XKZ  , NM,  NN  , JL6  , JM6  , JR  6 , JL  1 3 , JM  1 3 , 
5 JR  1 3 , JTT 1 , JTM 1 , JTT  2 , JTM2 , JTT  3 , JTM 3 , JTT 4 , JTM4 , JTT  5 , JTM  5 , JTT  7 , JTM7  , J 
6TT8, JTM8, JTT9, JTM9, JTT10, JTM10, JTT11 , JTM 1 1 , JTT 1 2 , JTM  1 2 , JTT  2 1 , JTM21 
7 , JTT 23 , JTM 23 , JTT 24 , JTM 2 4 , JTT  25 , JTM 2 5 

N-NN 

xco-xx 

XC 1 “X (N ) -2 . 0*XH 
XC  2*X (N )-XH 
XC3-X(N) 

ETA1-DIFF7 (XCO, XC 1 , ELL , CAY) 

ETA2-DIFF7 (XCO, XC 2 , ELL , CAY) 

ETA3**DIFF7  (XCO,  XC  3, ELL, CAY) 

CALL  PARAB  ( AAO , AA 1 , AA 2 ) 

RETURN 

END 


FUNCTION  3INT1  (A) 


A COMMONLY  USED  DEFINITE  INTEGRAL  APPEARING  IN  SEVERAL  EQ'NS. 
IS  EVALUATED  BY  THIS  ROUTINE.  TIIE  INTEGRAND  TAKES  THE  FORM  — 
P(X)*((ELL-X)/X)**.5.  — SEE  TRITE-UP. 


HEEDS  COMMON 
REAL  KA Y ( 5 0 ) 

COMMON  P (101  ) ,C(1  01  ) ,Q(101  ) ,11(101)  ,X(101  ) ,XSQ(101  ) ,FX(101  ),  ETN  1 (10 
1 1 ) , ETN 2 (101 ) ,ETP1 (101 ) ,ETP2( 101 ) ,YYU 1 ( 101 ) , YYU2 (101 ) ,TAU(101 ) ,ALP( 
2101 ) ,DEE ( 10  L ) , AA(50 ) ,KAY ,CAY ,DUM , A 1 , A 2 , AC A P , AL PH A , AMD A , ASQ, BB , BINT 

3 , CO  , CL  , CM,  DELTA  , DINT  , ECAP  , ELL  , EM,  E PS  ILN  , ET  A 1 , ET  A 2 , ET  A 3 , TE  F. , ON  E , T J 0 

4 , THR,  THETA,  TSA  VE  , XCO  , XC  1 , XC  2 , XC  3 , XH  , XKZ  , NM , NN  , JL  6 , JM  6 , JR  6 , JL  1 3 , JM  1 
53, JR13, JTT1, JTM1, JTT 2, JTM2, JTT 3, JTM3, JTT 4, JTM4, JTT 5, JTM5, JTT 7, JTM 7 
6, JTT 8, JTM8, JTT 9 , JTM9 , JTT l 0 , JTM 1 0 , JTT 1 l , JTM1 1 , JTT 12 , JTM  12 , JTT 21 , JTM 
72 1 , JTT 23 , JTM 23, JTT 24 , JTM  2 4 , JTT 2 5, JTM 2 5 

COMMON  / BLK  1 / V l BAR  , V 1 B AP.  1 , EC  AP  1 , YC  1 2 


ELL “A*  * 2 + 1 . 0 

XCO-O  . 0 
XC1-X(1  ) 
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XC2-X(2) 

XC  3-X (3 ) 

ETA1-P(1 )/SQRT(ELL-XCl) 

ETA2-P (2 )/SQRT (ELL-XC2) 

ETA3-P(3)/SQRT (ELL-XC3) 

CALL  PARAB  ( AAO , AA 1 , AA2 ) 

XC  l-X (NN-2 ) 

XC2-X(NN-1 ) 

XC  3-X (NN  ) 

ETA  1 -1 . O/SQRT ( ( ELL -XC 1 ) *XC  1 ) 

ETA2-1 . O/SQRT ( ( ELL-XC 2 ) *XC 2 ) 

ETA3-1 . O/SQRT ( ( ELL -XC 3 ) *XC 3 ) 

CALL  PARAB  ( BBO , BB  I , BB 2 ) 

NM-NN-4 

J-3 

DO  IOO  1-1, NM 

FX(I)“P(J)/SQRT((ELL-X(J))*X(J)) 

J-J  + l 

100  CONTINUE 

CALL  NTGRTE  (TCI) 

NM-3 

J-NN-2 

DO  110  I - 1 , NM 

FX(I )=C (J ) /SQRT ( (ELL-X(J ) )*X(J) ) 

J-J  + l 

110  CONTINUE 

CALL  NTGRTE  (TG2) 

BINT  1 -SIX (X (1 ) ,X(3 ) , AAO, AA 1 , AA2 )+TG 1+TG2+BB* FOUR (XC 1,BB0,BB1,3B2) 

RETURN 

END 


FUNCTION  BTG8  (DX) 
BTG8 


THE  DEFINITE  INTEGRAL  OF  P (X ) *LOG (DX+X** . 5 ) , WHICH  OCCURS  IN 
SEVERAL  EQUATIONS,  IS  EVALUATED  3Y  THIS  ROUTINE. 


MEEDS  COMMON 
REAL  KAY ( 50 ) 

COMMON  P(101),C(101),Q(101),H(101),X(101),XSQ(101),FX(101),ETN1(10 
1 1 ) ,ETN2  (101  ) , ETP  1 ( 10  1 ) , F.TP2(  101  ) , YYU  1 ( 101  ) , YYU2  (101  ) ,TAU  ( 101  ) , ALP( 
2101 ) ,DEE(10I ) ,AA(50) , KAY, CAY, A, Al, \2 , ACAP, ALPHA , AMD A , ASQ , BB , BINT , C 
30,  CL,  CM,  DELTA,  DINT  , ECAP,  ELL  , EM,  EPS  I LN  , ET  A 1 , ET  A 2 , ET  A 3 , TE  E , ON  E , T.J  0 , T 
4MR  , THETA  , TSAVE  , XC  0 , XC  1 , XC  2 , XC  3 , XI! , XKZ  , NM  , NN  , JL  6 , JM  6 , JR  6 , JL  1 3 , JM  1 3 , 
5 JR  1 3 , JTT  1 , JTM 1 , JTT2 , JTM2 , JTT  3 , JTM  3 , JTT4 , JTM4 , JTT  5 , JTM5 , JTT  7 , JTM  7 , J 
6TT8, JTM3, JTT9, JTM9, JTT  10, JTM10, JTT 1 1 , JTU1 L , JTT 12 , JTM 12 , JTT 21 , JTM 2 1 
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7, JTT23, JTM23, JTT24, JTM24, JTT25, JTM25 
COMMON  /BLK1/V1BAR, VI BAR  I , ECAP1 , YC12 
C 

IF  (DX)  100,100,110 
100  BTG8-0. 5*BTG5 (DX) 

GO  TO  140 
110  NM-NN-2 

DO  120  1-1 ,NM 

120  FX(I)“P(I)*CKLOG(D  X+X (I ) * * 0 . 5) 

CALL  NTGRTE  (TGI) 

J-NN-2 

NM-3 

DO  130  1*1, NM 

FX(I )-G (J )*CKLOG (D  X+X ( J ) * *0 . 5) 

J-J  + l 

130  CONTINUE 

CALL  NTGRTE  (TG  2 ) 

XCO-O.O 
XC  1 *=X  (NN  - 2 ) 

XC  2=X (NN- 1 ) 

XC 3=X (NN  ) 

ET A 1 “C  KLOG (D  X+XC 1**0. 5) 

ETA2-CKLOC (DX+XC2**0. 5) 

ETA3=CKLOG (DX+XC  3**0 . 5) 

CALL  PARAR  ( AAO , AA 1 , AA 2 ) 

BTG8=TG1+TG2+BB*F0UR (X(NN-2 ) , AA 0 , AA 1 , AA 2 ) 
140  RETURN 
END 


FUNCTION  BTG9  (A) 

C BTG9 

C 
C 
C 

C EVALUATES  THE  DEFINITE  INTECRAL  OF  P(X)/X**.5  OVER(O.l). 

C 

c 

REAL  KAY  (50) 

COMMON  P(101),G(101),Q(101),H(101),X(101),XSQ(101),FX (101), ETN 1(10 
1 1 ) ,ETN2 (101 ) ,ETP1 (101 ) ,ETP2( 101 ) ,YYU 1 (101 ) , YYU2(101 ) ,TAU ( 101 ) , ALP( 
2101 ) ,DEE(101 ) ,AA(50) , KA Y , C AY , DUMA , A 1 , A 2 , AC A P , ALPH A , AMD A , ASQ , 5R , B IN 
3T  , CO,  CL,  CM,  DELTA,  DINT  , F.CAP,  ELL  , EM , E PS  ILN  , ET  A 1 , ET  A 2 , ET  A 3 , TE  E , ON  E , TI 
40  , T1IR,  THETA  , TS  A VE  , XC  0 , XC  1 , XC  2 , XC  3 , Xll , XKZ  , NM , NN  , JL  6 , Jf  1 6 , JR  6 , JL  1 3 , JM 
5 13, JR  13, JTTl, JTM1, JTT2, JTM2, JTT3, JTM3, JTT4, JTM4, JTT5, JTM5, JTT  7, JTM 
6 7, JTT8, JTM8, JTT9, JTM9, JTT10, JTM10, JTT11 , JTM 11 , JTT  1 2 , JTM  1 2 , JTT 2 1 , JT 
7M21.JTT23, JTM 23, JTT24, JTM 24,  JTT25, JTM 2 5 
COMMON  /BLK1/V13AR , VI  BARI , ECAPl  , YC  12 
XCO-O. 0 
XC1-X(1 ) 

XC  2“X ( 2 ) 

XC3-X(3) 
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ETAI-P(1 ) 

ETA2-P (2 ) 

ETA3=P(3) 

CALL  PARAB  ( AAO , AA 1 , AA 2 ) 

NM-NN-4 

J-3 

DO  100  1=1, NM 

FX(I)-P(J)/SQRT(X(J)) 

J-J+l 

100  CONTINUE 

CALL  NTCRTE  (TGI) 

NM-3 

J-NN-2 

DO  110  1=1,3 

FX(I)=C(J)/ SQRT (X ( J ) ) 

J=J  + 1 

110  CONTINUE 

CALL  NTCRTE  (TC2 ) 

BB 1=SQRT (X  (3 ) ) 

BB2=SQRT (X (NN-2 ) ) 

BTC9=SIX(X( 1 ) , X ( 3 ) , AAO , AA 1 , A A 2 ) +T  G 1 +T  G 2+B  B* ( AT AN (BB1/3B2)-BB1*BB2) 

RETURN 

END 


SUBROUTINE  NTCRTE  (TG) 
NTGR 


THIS  SPECIAL  INTEGRATION  ROUTINE  IS  USED  REPEATEDLY  BY  THE 
PROCRAM.  HIGH-ORDER  NET  TON-COTES  FORMULAS  ARE  USED 


REAL  KAY (50) 

COMMON  P(101)  ,G(  101  ) ,Q(10i  ) ,11(101  ) ,X(  101  ) ,XSQ(10l  ) ,FX(  101  ) , ETN  1(10 
1 1 ) ,ETN2 (101 ) ,ETP1  (101  ) ,ETP2 ( 101 ) , YYU 1 ( 101 ) , YYU2 ( 101 ) ,TAU (101 ) ,ALP ( 
21  01  ) ,DEE  (101  ) , AA(50)  , KAY , CAY  , A , A 1 , A2  , ACAP  , ALPHA  , A'IDA  , ASQ  , BB  , BINT  , C 
30 , CL , CM, DELTA , DINT , ECAP , ELL , EM, EPSILM , ETA  1 , ET A 2 , ET A 3 , TE E , TC  1 , TG 2 , T 
4G 3, THETA , TS AVE , XCO, XC 1 , XC  2 , XC  3 , XH  , XKZ , NM, NN , JL 6 , JM6 , JR 6 , JL 1 3 , JM1 3 , 
5 JR  1 3 , JTT 1 , JTM 1 , JTT  2 , JTM2, JTT  3 , JTM 3 , JTT  4 , JTM4 , JTT  5 , JTM  5 , JTT  7 , JTM  7 , J 
6TT8, JTM8, JTT9, JTM9, JTT10, JTM10, JTT11 , JTM1 L , JTT 1 2 , JTM 1 2 , JTT 2 1 , JTM 21 
7 , JTT 23 , JTM 23 , JTT 24 , JTM 24 , JTT 2 5 , JTM 2 5 

DIMENSION  YNN(IO) 

TC  = 0 . 0 

IF  (NM-1)  100,100,120 
100  I RITE  (6,110) 

110  FORMAT  (54(10  ***ERROR  RETURN  FROM  NTCRTE — NO.  OF  CRIO-PTS.  L.T.l) 

GO  TO  390 
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IF  (NM-8 ) 140,130,140 

H7-1 

LN-0 

KN-0 

CO  TO  230 

ND-NM-1 

N7-0 

IF  (ND-12)  170,160,160 
ND-ND-7 
N7-N7  + 1 
GO  TO  150 

IF  (ND-ll)  180,220,180 
IF  (ND-7 ) 190,200,210 
KN-ND 
LN-0 

GO  TO  230 
N 7-N  7+1 
KN  =0 
LN-0 

CO  TO  230 
KN  “4 

LN-ND-KN 
GO  TO  230 
KN-6 
LN  =5 
N STRT  = 1 
NN7-N7 

IF  (NN  7 ) 250,280,250 
IS/ 1-7 
J -N  STRT 
DO  260  1=1,8 

YNN (I )=FX(J ) 

J-J  + l 
CONTINUE 

ASSICN  270  TO  KLOC 
GO  TO  400 
TG-TC+TG1NT 
NSTRT-N’STRT  + 7 
NN7-NM  7-1 
GO  TO  240 

IF  (KN)  290,320,290 

IS J I =KN 

KK-KN+1 

J-NSTRT 

DO  300  1=1, KK 

YNN  (I ) = FX(J ) 

J-J  + l 

COKT INUE 

NSTRT-NSTRT+KN 

ASSICN  310  TO  KLOC 

GO  TO  400 

TC-TC+TCINT 

IF  (LN)  330,360,330 

I S' I I =LN 

KK-LN+1 

J-NSTRT 


250, 280, 250 


KLOC 


1 1 , KK 


o no  noon 


* ..  — A 
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DO  340  I-l.KK 
YNN(I)-FXCJ) 

J-J  + I 

340  CONTINUE 

NSTRT-NSTRT+LN 
ASSIGN  350  TO  KLOC 
GO  TO  400 
350  TC-TC+TGINT 
360  IF  (NM-NSTRT ) 370,390,370 
370  M RITE  (6,380) 

380  FORMAT  (67H0  ****ERROR  IN  FINAL  VALUE  OF  STEP  COUNTER  (NSTRT) — NT 
1GRTE  ROUTINE) 

GO  TO  390 
390  RETURN 

400  IF  ( I SM  I ) 410,410,430 
410  '1  RITE  (6,420) 

420  FORMAT  (47H0  *****ERR0R  IN  ISM  SETTING  IN  NTCRTE  ROUTINE) 

GO  TO  390 

430  IF  (ISM  1-2 ) 470,480,440 
440  IF  (ISM  1-4)  490,500,450 
450  IF  (ISM  1-6  ) 510,520,460 
460  IF  (ISM  1-7  ) 370,530,370 
4 70  TGINT= (YNN ( 1 )+YNN (2 ) ) *XH/2. 0 
GO  TO  540 

4 80  TGINT=  (YNN  ( 1 ) + 4. 0*YNH  (2  )+YNN  (3  ) ) *X!I/3. 0 
GO  TO  540 

4 90  TGINT= (YNN ( 1 ) + 3. 0 * Y N N (2 ) + 3. 0*YNN (3 )+YNN (4 ) ) *3. 0*Xll/8  . 0 

GO  TO  540 

5 00  TGINT  = ( 7 . 0 *YNN ( 1 ) + 32 . 0*YNN ( 2 ) + 1 2 . 0 * Y N N (3 ) + 32 . 0*YNN (4 )+7 . 0*YNN (5 ) ) * 

14. 0 *XH / 90 . 0 
GO  TO  540 

510  TGINT  = (19. 0 *Y N N (1 ) + 75. 0 * Y N N (2)  + 50. 0*YNU (3)  + 50. 0*YNN (4 ) + 7 5. 0*YNN (5) 
1+19. 0*YNN (6) ) *5 . 0*XH/288. 0 
GO  TO  540 

520  TGI NT = (4  1 . 0*YNN  ( 1 ) + 216. 0*YNN  (2)  1-2  7.  0 * Y N N ( 3 ) + 2 7 2 . 0 *Y  NN  ( 4 ) + 2 7 . 0*YNN  ( 
15 )+2 16. 0*YNN (6 )+4 1 . 0*YNN(7))*6. 0*XH/840. 0 
GO  TO  540 

5 30  TGI  NT  = (75 1 . 0*YNN ( 1 ) + 35  7 7. 0 * Y N N (2 ) + l 323. 0*YNN ( 3 ) + 2 989 . 0*YNN(4)+29S9 
1 . 0*YNN (5)  + l 323. 0*YNN (6 ) + 3577 . 0*YNN (7 ) + 7 51 . 0*YNN (8) ) *7. 0*XH/  1 7280. 0 
540  CONTINUE 

GO  TO  KLOC,  (270,310,350) 

END 


FUNCTION  BTG4  (X l , X2 , BB 0 , RB 1 , BB 2 ) 
BTG4 


THIS  SUBSIDIARY  ROUTINE  EVALUATES  THE  DEFINITE  INTEGRAL  OVER 
(0,211)  WHOSE  INTEGRAND  IS  (BH  1 *X  + BB2*X**2  ) *LOC  (X  ) — REQUIRED 


... 
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C BY  THE  ROUTINE  BTG5. 

C 

C 

c 

C NO  COMMON. 

C 

FLOW  -0. 0 

IF  (X1.NE.0.0)  FLOW =BB1*X1 **2*(AL0G(X1 )-0. 5)/2. 0+BB2*Xl **3*(ALOG(X 
11)  — 1. 0/3. 0)/3. 0+BB  0*X1 * ( ALOG (Xl)-l.O) 

BTG4-BB  l*X2**2*(AL0G(X2)-0. 5 ) / 2 . 0+BB2*X2**3*(ALOG(X2 )-l . 0/3. 0 ) / 3 - 0 
1+3B 0*X2  * ( ALOC  ( X 2 ) — 1 . 0)-FLO.J 
RETURN 
END 


FUNCTION  BTC  5 (A) 

C 

C 

C 

C THE  DEFINITE  INTEGRAL  OVER  (0,1)  WITH  INTECRAND  P(X)LOC(X) 

C IS  EVALUATED. 

C 

c 

c 

c 

C NEEDS  COMMON. 

C 

REAL  KAY (50) 

COMMON  P(101)  , G ( 1 0 1 ) ,0(101  ) ,11  (101  ) , X ( 1 0 l ) , XS  Q ( 1 0 1 ) , FX  ( 1 0 1 ) ,ETN1  (10 
11 ) ,ETN2(101 ) ,ETP1 (101 ) ,ETP2(101 ) , YYU1 (101 ) , YYU2  (101 ) ,TAU(101 ) ,ALP( 
2101 ) , DEE (101 ) ,AA(50) , KAY , CAY , DUMA , A 1 , A 2 , AC AP , ALPHA , AMD A , ASQ , B B , B IN 
3T,  CO,  CL,  CM,  DELTA,  DINT,  ECAP,  ELL  , EM , EPS  ILN  , ETA  1 , ETA2  , ETA  3 , TE  E , ON  E , T-J 
40, THR, THETA , TSAVE , XCO, XC 1 , XC2, XC3 , XH , XKZ , NM.NN , JL 6 , JM6 , JR 6 , JL  1 3 , JM 
51 3, JR  13 , JTT 1 , JTM1 , JTT  2 , JTM  2 , JTT  3 , JTM  3 , JTT  4 , JTM4 , JTT  5 , JTM 5 , JTT  7 , JTM 
6 7 , JTT  8 , JTM8 , JTT  9 , JTM 9 , JTT 10 , JTM 10 , JTT 1 1 , JTM 1 1 , JTT 1 2 , JTM 1 2 , JTT  2 l , JT 
7M21 , JTT  2 3 , JTM 2 3 , JTT 24 , JTM 24 , JTT 2 5 , JTM 2 5 
COMMON  /BLK1 /VI BAR , VlBARl , ECAP1  , YC12 
C 

NM=NN -4 
J = 3 

DO  100  1=1, NM 

FX(I )=P (J )*CKLOG (X (J ) ) 

J=J+1 

100  CONTINUE 

CALL  NTCRTE  (TGI) 

XC0=0. 0 
XC1=X(1  ) 

XC  2=X ( 2 ) 

XC  3 = X ( 3 ) 

ETA1=P(1  ) 

ETA2-P  (2  ) 

ETA3-P  (3 ) 


f-wov. .y.. 
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CALL  PARAB  ( AAO , AA 1 , AA  2 ) 

NM-3 

J-NN-2 

DO  110  I-l.NM 

FX(I )«G (J ) *CKLOG (X  (J  ) ) 

J-J  + l 

110  CONTINUE 

CALL  NT CRT E (TG  2 ) 

BTG5-BTG4(X(1 ) ,X(3) , AA 0 , AA 1 , AA 2 ) +TG 1 +TG 2+B B * BTG6 (XH ) 

RETURN 

END 


FUNCTION  BTG6  (XH) 


C 

C 

C THE  DEFINITE  INTEGRAL  OVER  (1-2H,1)  WITH  INTEGRAND 

C L0G(X)*( 1-X)**. 5 

C IS  EVALUATED.  - REQUIRED  BY  BTG 5 . 

C 

c 

C 

C 

c 

C NO  COMMON. 

C 

XI ■ ( 2 . 0*XH  ) **  1 . 5 
X2  *»  ( 2 . 0*XH)**0. 5 
X3»(X1-1. 0 ) *CKLUG (X2-1 . 0) 

X3=X3+ (X 1+1 . 0)*CKLOG(X2+1 . 0) 

X3-X3-2. 0*X1  /3. 0 
X3-X3-2. 0*X2 
BTG6=2 . 0*X3/3. 0 
RETURN 
END 


SUBROUTINE  PARAB  ( CON  0 , CON  1 , CON  2 ) 

C PARB 

C 
C 
C 

C PUTS  A QUADRATIC  OF  THE  FORM 

C CON  O+CON 1 * ( X-XC  0 ) +C0N  2 * ( X-XC 0 ) **  2 

C THROUGH  THE  POINTS  (XC 1 , ET A 1 ) , ( XC 2 , ET A 2 ) , ( XC 3 , ET A 3 ) . 

C 


u o 
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REAL  KAY (50) 

COMMON  P(I01),G(101),Q(101),H(10I),X(10I),XS Q(IOI), FX (101)  , ETN  1(10 
11 ) ,ETN2 (101 ) , ETP1 (101 ) ,ETP2 (101 ) , YYU 1 (101 ) , YYU2( 101 ) ,TAU (101 ) , ALP( 
2101), DEE  (101),  A A (50),  KAY  , CA  Y , A , A 1 , A2  , ACA  P , AL  PH  A , AMDA  , AS<) , BR  , BItIT  , C 
30 , CL , CM ,DELTA,DINT,ECAP,ELL,EM,EPSILN,ETAl,ETA2,ETA3,TEE,ONE,TJO,T 
4HR  , THETA  , TSAVE  , XC  0 , XC  1 , XC  2 , XC  3 , XH  , XKZ  , NM , NN  , JL  6 , JM6  , JR  6 , JL  1 3 , JU 1 3 , 
5 JR  13, JTT1, JTM1, JTT2, JTU2, JTT3, JTM3, JTT  4, JTM4, JTT5, JTM5, JTT7, JTM7,  J 
6TT8, JTM8, JTT9, JTM9, JTT 10, JTM10, JTT 1 1 , JTM11 , JTT 12 , JTM12 , JTT21 , JTM21 
7, JTT23, JTM23, JTT24, JTM24, JIT  25, JTM25 

DENM 1 - (XC  2-XC 1 ) * (XC  3-XC 1 ) 

DENM2-(XC 1-XC2)* (XC3-XC2) 

DENM 3” (XC 1 -XC  3 ) * (XC 2-XC  3 ) 

CON  0-ETA 1 * ( XC  2-XC  0 ) * ( XC  3-XC  0 ) /DENM 1+ETA 2 * ( XC 1 -XC  0 ) * (XC  3-XC 0 ) /D ENH 2 
1+ETA 3* (XC 1-XC0)*( XC  2-XCO ) /DENM 3 

CON  1 -ETA  1* (2. 0*XC0-XC2-XC3) /DENM 1+ETA2* ( 2 . 0 *XC 0-XC  1 -XC 3 ) /DENM2+ETA 
1 3* ( 2 . 0*XC0-XC 1 -XC  2 ) /DENM 3 
CON2-ETA1 /DENM1+ETA2/DENM2+ETA3/DENM3 
RETURN 
END 


FUNCTION  SIX  (XI , X2 , SFO , SF 1 , SF 2 ) 

C EVALUATES  THE  INTEGRAL  FROM  0 TO  H OF 
C SFO+SF 1*T+SF2*T**2 ) /SORT (T ) . 

C 

C 

C 

C 

SIX=0. 06G666667*(30. 0*SF 0+X2 * ( 1 0 . 0 *S  F 1 +6  . 0 *S  F 2 *X 2 ) ) *SQRT  (X2  )-0. 06  6 
1666667* (30. 0*SF0+X1 * ( 10. 0*SF 1+6. 0*SF2*X1 ) ) *S0RT (XI ) 

RETURN 

END 
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RETURN 

110  CKLOG-ALOG (-CARGLE) 
RETURN 

120  CKLOG-ALOG (GARGLE) 
RETURN 
END 


FUNCTION  FOUR  (H , CFO , CF 1 , CF 2 ) 

FOUR 

EVALUATES  THE  EXPRESSION  OBTAINED  FOR  THE  INTEGRAL  FROM  H TO  1 OF 
(CFO+CF 1*T+CF2*T**2 )*SQRT ( 1-T ) . 


FOUR-O  . 0190  47619*05. 0*CF0  + 7 . 0*CF  1*(3. 0*H+2. 0 ) +C  F 2 * ( H * ( 1 5 . 0*11  + 12 . 0 
1 ) + 8. 0) ) *(1 . 0-H ) * * 1 . 5 
RETURN 
END 


SUBROUTINE  REGUL  ( FUNNY  , HERE  , GAUCHE  , DROIT  , SMALL  , ROOT  , PCON  , NS  I Cll) 
REGU 


FINDS  THE  ROOT  OF  THE  EQUATION  — 

FUNNY  (X,  PCON  ,NSJ  Cll) -HERE-0 
BY  THE  METHOD  OF  REG1JLUS  FALSI. 


FUNNY  MUST  APPEAR  ON  AN  F-CARD  IN  THE  CALLING  ROUTINE 
ASSUMES  FUNNY  (X  , PCON  , IIS  I CH  ) 

ICNT-0 

VR-FUNNY (DROIT , PCON , NS  I CH) 

VL-FUNNY (GAUCHE, PCOM.NSJ CH) 

IF  ( ( VL-HERF,  ) * ( VR-HERE  ) ) 100,240,260 
100  EMN-DROIT 

110  IF  (IC NT-50)  120,160,160 

120  EM- (DROIT -GAUCHE  )*  (HERE-VL  ) / ( VP.-VL  )+CAUCHE 
130  IF  ( ABS (EM-EMN )- SMALL ) 140,170,170 
140  ROOT -EM 
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GO  TO  110 
210  GAUCHE-EM 


VL-VM 
GO  TO  200 

220  IF  ( VM-HERE ) 230,140,230 

230  ROOT  =GAUC1IE 
GO  TO  150 

240  IF  (VL-HERE)  250,230,250 

250  ROOT-DROIT 
GO  TO  150 

260  '.1  RITE  (6,  270) 

270  FORMAT  (64II0  *****ERROR  EXIT  FROM  REGUL  — END-PTS . INPUT  DO  NOT  BRA 
1CKET  ROOT) 

GO  TO  150 

280  WRITE  (6,290)  GAUCHE  , EM  , DROIT  , VL  , VR  , HERE 

290  FORMAT  (37HO  ****  REGUL  ITERATION  COUNT  PAST  100/6X,  14HLEFT-HAND  P 
IT. -El  4. 6, 6X,  13HLAST  APPROX. =E  1 4 . 6 , 6X , 6X , 13HRT . -HAND  PT . =E 1 4 . 6 / 6X , 1 
25HFUNCT. -AT  LEFT  = E14. 6,6X,  14HFUNCT.  AT  RT.-E14.6.6X,  18HINPUT  FUNCT 
3.  VAL . =E 1 4 . 6 ) 

STOP 

END 


C 

C 

C 

C 

C 

C 

C 

C 

C 

C 

C 

C 


FUNCTION  BTG1  (A) 


A COMMONLY  USED  DEFINITE  INTEGRAL  APPEARING  IN  SEVERAL  EQ'NS. 
IS  EVALUATED  BY  THIS  ROUTINE.  THE  INTEGRAND  TAKES  THE  FORM  — 

P (X ) * ( ( ELL-X ) /X ) * * . 5 . — SEE  WRITE-UP. 


NEEDS  COMMON 
REAL  KAY (50) 

COMMON  P ( 1 0 1 ) , G ( 1 0 1 ) , Q ( 1 0 1 ) , H ( 1 0 1 ) , X ( 1 0 1 ) , XS  Q ( 1 0 1 ) , FX ( 1 0 1 ) , ET  N 1 ( 1 0 
1 1 ) ,ETN2(101 ) ,ETP1 (101) ,ETP2(101 ) , YYU  1 ( 1 0 1 ) , Y Y U 2 ( 1 0 1 ) ,TAU(101) ,ALP( 

2 l 01  ) , DEE  ( 10  l ) , AA(  30  ) , KAY  , CAY  , DU ‘I , A1  , A 2 , AC.AP,  ALPHA,  AMD  A , ASy , BB , BINT 

3 , CO  , CL  , CM  , DELT  A , DI  NT  , EC  A P , ELL,  F.M  , E PS  I LN  , ET  A 1 , ET  A 2 , ET  A 3 , TE  E , ON  E , Tf  0 

4 , TUR  , THETA  , TS  A VE  . XCO  , SC  l , XC  2 , XC  3 , XH  , XKZ  , KM  , NN  , JL  6 , JM 6 , JR  6 , JL  1 3 , JM  1 


o o 


r 


wpj  ' 


i 
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5 3 , JR  1 3 , JTT  1 , JTM  1 , JTT 2 , JTM2  , JTT  3 , JT113  , JTT  4 , JIM 4 , JTT  5 , JIM 5 , JTT  7 , JTM7 
6, JTT8, JTM8, JTT9, JTM9, JTT 10, JTM10, JTT 1 1 , JTM11 , JTT  1 2 , JTM 1 2 , JTT  2 1 , JIM 
721, JTT 2 3, JIM 2 3, JTT 24, JTM 24, JTT 25, JTM 2 5 
COMMON  / BLK 1 / V 1 3 AR , VI BAR  1 , ECAPl , YC  12 


ELL-A**2+1 . 0 
XCO-O. 0 
XC1-X(1  ) 

XC2-X(2  ) 

XC 3-X (3 ) 

ETAl-P (1 )*SQRT (ELL-XC  1 ) 

ETA2-P (2 )*SQRT (ELL-XC 2) 

ETA3-P (3 )*SQRT (ELL-X(3  ) ) 

CALL  PAP.AB  ( AAO  , A A 1 , A A 2 ) 

XC 1 »X (NN  - 2 ) 

XC  2-X (NN - 1 ) 

XC  3-X (NN  ) 

ETA1-SQRT ( (ELL-XC 1) /XC  1 ) 

ETA2-SQRT ( (ELL-XC 2) /XC2) 

ETA3-SQRT  ( ( ELL-XC  3 ) /XC  3 ) 

CALL  PARAB  ( BBO , BB 1 , BB 2 ) 

NM-NN-4 

J-3 

DO  100  I-l.NM 

FX(I)=P(J)*SQRT((ELL-X(J))/X(J)) 

J-J  + l 

100  CONTINUE 

CALL  NTGRTE  (TGI) 

NM-3 

J-NN-2 

DO  110  I-l.NM 

FX ( I ) =C (J )*SQRT ( (ELL-X(J ) ) /X(J  ) ) 

J-J+l 

110  CONTINUE 

CALL  NTCRTE  (TG  2 ) 

BTGl»SIX(X(l),X(3) , AAO, AA1 , AA2)+TG  1+TG 2 + BB* FOUR (XC  1,BB0,BB1,BB2) 

RETURN 

END 


C 

C 

C 

c 

c 

c 

c 


FUNCTION 

FTJX 


FT  JX  (II, XX) 


THIS  SUBSIDIARY  ROUTINE  EVALUATES  THE  DEFINITE  INTEGRAL  JHOSE 
INTEGRAND  IS  P ( U ) * LOC ( U -XX ) FOR  0 LT  XX  LE  1 — U IS  THE  DUMMY 
VARIABLE  OF  INTEGRATION. 


1^  I 
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C NEEDS  COMMON 

REAL  KAY (50 ) 

COMMON  P(10i  ) ,G(  101  ) ,Q(I01)  ,11(101  ) ,X(10i  ) ,XSQ(  101  ) ,FX(  101  ) ,ETN1  (10 
II  ) ,ETM2(  101  ) ,ETP1  ( 101  ) ,ETP2(  101  ) , YYIJ  1(101)  , YYU2  (101  ) ,TAU  (101  ) , ALP( 
2101)  ,DEE(101)  , AA  (50),  KAY,  CAY,  A,  Al  , A2  , ACAP  , ALPHA  , AMD  A , ASl)  , BB  , 3INT  , C 
30 , CL, CM, DELTA, DINT , ECAP, ELL , EM , E PS ILN , ET A 1 , ET A2 , ET A 3 , TEE , ON E , TW  0 , T 
4HR  , THETA  , TS A VE  , XC  0 , XC  1 , XC  2 , XC  3 , XH  , XKZ  , NM , NN  , JL  6 , JM 6 , JR 6 , JL  1 3 , JM  1 3 , 
5 JR  1 3 , JTT 1 , JTM 1 , JTT  2 , JTM2 , JTT  3 , JTM  3 , JTT4 , JTM4 , JTT  5 , JTM5 , JTT  7 , JTM7 , J 
6TT8, JTM8, JTT9, JTM9, JTT10, JTM10, JTT11 , JTM 11 , JTT  l 2 , JTM  1 2 , JTT 2 l , JTM21 
7, JTT 23, JT  123, JTT 24, JTM 24, JTT 25, JTM2 5 
COMMON  /BLK1/V13AR, V13AR1 , EC API  , YC12 
IF  ( I I- 1 ) 180,130,100 
100  IF  ( 1 1-2  ) 190  , 190  , 1 10 
1 10  IF  ( 1 1 -3 ) 200,200,  120 
120  IF  (II-NN)  130,260,180 
130  IF  ( 1 1 - (NN-4 ) ) 210,  220,  140 
140  IF  (II— (NN-3) ) 230,230,150 
150  IF  (II-(NN-2) ) 240,240,160 
160  IF  (II-(NN-l))  250,250,170 
170  WRITE  (6,280) 

STOP 

180  FTJX-EJ1 (1 ,NN , XX , P) 

RETURN 

190  XC  0=XH+X ( 1 ) 

XC1«X(1) 

XC2«X(2) 

XC  3=X ( 3 ) 

ETA  1 *=P  ( 1 ) 

ETA  2=P ( 2 ) 

ETA3=P (3) 

CALL  PAR A3  ( AAO , AA  1 , AA2 ) 

NN2-NM-2 

FT  JX  = E J 3 (AAO  , AA  2 , XH  ) + E J 1 ( 3 , NN  2 , XH  , P ) +E  J 1 (N  M 2 , NN  , XH  , G ) +B  3 * E J 4 ( XH  , XH 
1) 

GO  TO  270 

200  XC0=2. 0*XH+X(1 ) 

XC1=X(1  ) 

XC  2*=X  ( 2 ) 

XC  3=X ( 3 ) 

ET  A 1 =P  ( 1 ) 

ETA2-P (2 ) 

ETA  3 = P ( 3 ) 

CALL  PARAB  ( AAO , AA 1 , AA2 ) 

XC 1 =X ( 3 ) 

XC  2=X  ( 4 ) 

XC 3=X  (5  ) 

ETA  1 *=P  ( 3 ) 

ETA2-P(4) 

ETA  3“P ( 5 ) 

CALL  PARAB  ( BCO , BB  1 , BB 2 ) 

NN2-NN-2 

FTJX-EJ5 ( AAO, AA1 , AA 2 , BB 0 , BB 1 , 3B 2 , XH ) + 

1 ,C)  + BB*EJ4  (XII  , XX) 

GO  TO  270 
2 10  I1L-II-1 
I1R-I 1+1 
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XCO-XX 

XCl-X(IIL) 

XC2-XX 

XC3«X(IIR) 

ETAl-P(IIL) 

ETA2-P (II) 

ETA3-P  (HR) 

CALL  PARAB  ( AAO , AA I , AA2 ) 

NN2-NN-2 

FTJX-EJ1  (1 , IIL , XX ,P)  + EJ3( AAO , A A 2,XH)  + EJI ( I I R , NN 2 , XX , P )+E J 1 (NN 2 , NN , 
1XX, C)+33*EJ4 (XH.XX) 

GO  TO  270 
220  NN6-NN-6 
NN2-NN-2 
NN5-NN-5 
NN3-HN-3 
XCO-XX 
XC 1-X (NN  6 ) 

XC  2-X ( NN  5 ) 

XC3-XX 
ETA1-P ( H N 6 ) 

ETA2-P (NN5 ) 

ETA3-P (II) 

CALL  PARAB  ( AAO , AA 1 , AA 2 ) 

XC  1-XX 
XC  2-X (NN  3 ) 

XC  3-X (NN  2 ) 

ETA1-P (I  I ) 

ETA2-P (NN3) 

ETA3-P (NN  2 ) 

CALL  PARAB  ( BBO , BB 1 , BB 2 ) 

FTJX-EJ1  ( 1 , tiN  6 , XX  , P ) + E J 5 ( AAO  , AA 1 , AA  2 , BB  0 , BB 1 , BB  2 , XH ) +E J 1 (NN  2 , NN , XX 
1 , C)+33*E J4 (XH , XX) 

CO  TO  270 
230  XCO-XX 

XC1-X(II-1 ) 

XC2-XX 
XC  3-X (I I + I ) 

ETA1«P(II-1 ) 

ETA2-P (II ) 

ETA3-P (II+l ) 

CALL  PARAB  ( AAO , AA 1 , AA 2 ) 

IIL-II-1 

IIR-II+1 

FTJX-EJ1 (1 ,IIL,XX,P)+EJ3( AAO ,AA2,XH)+EJ1 ( I I R , NN , XX , G ) +B B* E J 4 (XH , XX 
1) 

GO  TO  270 
240  XCO-XX 

XC 1-X (1 1-2 ) 

XC  2-X ( I I - 1 ) 

XC3-XX 

ET  A 1 -P  ( 1 1 - 2 ) 

ETA2-P(II-1 ) 

ETA3-P (II ) 

CALL  PARAB  ( A A 0 , AA 1 , AA 2 ) 

XC  1-XX 
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XC2«X(II+1 ) 

XC3-X(II+2) 

ETA1-C (II ) 

ETA2-G(II+1 ) 

ETA3*G ( I 1+2 ) 

CALL  PARAB  ( B B 0 , B B 1 , BB 2 ) 

IIL-II-2 

FTJX-EJ1 (I , IIL.XX, P)+EJ5(AA0,AA1,AA2,BB0,BB1,BB2,XH)+BB*EJ4(XH,XX) 
GO  TO  270 
250  XCO-XX 

XC1-X(II-I ) 

XC2-XX 

XC3«X(NN) 

ETAl*C (NN-2 ) 

ETA2*>G (NN  — 1 ) 

ETA3-G  (Nil ) 

CALL  PARAB  ( AAO , AA 1 , AA2 ) 

IIL-II-l 

FTJX-EJ1  (1 , IIL , XX,P)+EJ3(AA0,AA2,XH)+BB*EJ4 (XH.XX) 

GO  TO  270 
260  XC0*X(NN) 

XC 1-X (UN-2 ) 

XC2«X(HN-1 ) 

XC  3“X  (NN  ) 

ETA1-G (NN-2  ) 

ETA2=G (NN- 1 ) 

ETA  3*C (NN) 

CALL  PARAB  ( AAO , AA 1 , AA2 ) 

NN  2=NN-2 

FT JX  = E J 1 (1 , NN2, XX, P )+E J5 (AAO , AA 1 , AA2,0. 0,0. 0,0. 0 , XH  ) +B  B * E J 4 ( XH  , XX  ) 
270  RETURN 

280  FORMAT  ( 1H0 , 10X, 9H *********  , 27H ERROR  IN  FOPTAB  — FTJX  ROUT.) 

END 


FUNCTION  F1FUNC  (XNU.A) 
FIFN 


SPECIAL  FUNCTION  NEEDED  FOR  THE  NOSE  SINGULARITY  CASE  (IN  THE 
CALCULATION  OF  THE  CAVITY)  AND  FOR  THE  OFF-DESICN  CALCULATION. 


DOUBLE  PRECISION  D ASQ , DMCS Q , DB ET 1 , DB ET 2 , DC AM  1 , DC AM  2 

AS  Q”A*  * 2 

DASQ-ASQ 

DMCSQ-A*DS  QRT  (DASQ+l . 0) 

DBET1—  DSQRT((DMCSQ-DASQ)/2. 0D0) 

DBET2--DSQKT ( (DMGSQ+DASQ) /2 . ODO) 


noon 


-61- 

HYDROFOIL  DESIGN  PROGRAM 


16  June  1977 
JF: jep 


C 


DGAM1- (A*DBET2+DBETl/2.0)/0MCSQ 

DC AM 2“ (A*DBETl-DBET2/2. 0) /DMGSQ 

BET1-DBET1 

BET2-DBET2 

GAM  1 -DC AMI 

GAM2-DCAM2 

DELT l-A*DBET2/2. 0 /DMGSQ 
DELT2-A*DBET 1/2. 0 /DMGSQ 
OMGSQ-DMCSQ 

F 1FUNC=XNU*SQRT  (XNU*  (XNU+1 . 0))  / (XNU**2+ASQ) 

F1FUNC-F1FUNC+3ET1*ALFUNC(XNU,A,0MCSQ,BET1,BET2,DELT1,DELT2,GAM1,G 
1 AM  2 ) / ( 4 . 0*OMCS  Q ) 

F IFUNC-F 1FUNC+BET  2*TEFUNC (XNU , A, DELT I , DELT2 , CAM  I , GAM  2 ) / ( 2 . 0*OMCSQ) 

F1FUMC-0. 5*F 1FUNC/ AS  Q 

RETURN 

END 


FUNCTION  EJ1  ( I L , IR , XX , F ) 

NEEDS  COMMON  — USES  NTCRTE 
REAL  KA  Y ( 5 0 ) 

COMMON  P(101  ) ,G(  101  ) ,Q(101  ) ,11(101  ) ,X(101  ) ,XSQ(  101  ) ,FX(10l)  ,ETN1  (10 
11  ) , ETN 2(101)  ,ETP1  (101 ) ,ETP2(101 ) , YYU 1 (101 ) , YYU2 (101 ) ,TAU( 101 ) ,ALP ( 
2101 ) ,DEE (101 ) , AA(50) , KAY , CAY , A , A 1 , A2 , ACAP , ALPHA , AMDA , ASQ , BB , BINT , C 
30,CL,CM,DELTA,DIHT,ECAP,ELL,EM,EPSILN,ETA1,ETA2,ETA3,TEE,ONE,TJO,T 
4HR  , THET  A , TSAVE  , XCO  , XC  1 , XC  2 , XC  3 , XII , XKZ  , NM  , NN  , JL  6 JM6  , JR6  , JL  1 3 , JM  1 3 , 
5 JR  13 , JTT 1 , JTM1 , JTT  2 , JTM 2 , JTT  3 , JTM3 , JTT  4 , JTM4 , JTT  5 , JTM5 , JTT 7 , JTM7 , J 
6TT8, JTH8, JTT9, JTM9, JTT10, JTM10, JTT11 , JTM 1 1 , JTT 1 2 , JTM 1 2 , JTT 2 1 , JTM 21 
7 , JTT 23 , JTM 23 , JTT  24 , JTM 24 , JTT 2 5 , JTM 2 5 
COMMON  /BLK1 /VI  BAR , VI  BARI , ECAP1 , YC 12 
DIMENSION  F(101) 

NM-IR-IL+1 

J-IL 

DO  100  I-l.NM 

FX(I )-F( J)*CKLOG (X(J)-XX) 

J-J  + l 

100  CONTINUE 

CALL  NTCRTE  (EJ1) 

RETURN 

END 


FUNCTION  E J 3 (AZ.A2.XH) 

EVALUATES  THE  EXPRESSION  FOR  THE  DEFINITE  INTEGRAL  OF 
P (T  ) *LOG  (T-XX  ) FROM  XX-H  TO  XX+11. 


V u 
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(3-PT.  INTEGRAL  ACROSS  LOG.  SINGULARITY.  ) 

EJ3-2. 0*AZ*Xll*(CKLOC  (XH)-1 . 0)  + 0.  66 66 6 6 6 7 * A2 *XH ** 3 * ( CKLOC  (XH)-O.  333 
133333) 

RETURN 

END 


FUNCTION  EJ4  (XH,XX) 

C 

C EVALUATES  THE  DEFINITE  INTEGRAL  OF  ( 1-T ) ** . 5*LOG (T-XX ) OVER 

C (1-2H,1)  AS  OBTAINED  IN  CLOSED  FORM  (T  THE  VARIABLE  OF  INTEGRATION 

C 

Q 

C 

C INTEGRATE  SQ.  ROOT  TERM  FROM  1-2H  TO  1. 

C 

IF  ( AB  S(XX-(1 . 0-2. 0*XH) )-5. 0E-07)  120,  120,  100 
100  EJ4=(  (2. 828428*XH**1. 5-(l . 0-XX)**l . 5)*CKLOG ( 1 . 4142 14*SQRT (XH)-SQRT 
1 (1. 0-XX) ) + (2. 82  8428*XH**1  . 5+(l . 0-XX)**l . 5)*CKLOG(l . 4142 14*SQRT  (XH) 
2+SQRT (1 . 0-XX) )-l. 88561 9*XH** 1 . 5-2. 828428*SQRT (XH)*( 1 . 0-XX) ) *0.6666 
36667 

110  RETURN 

120  EJ4=(  (2. 82842  8 *XH  * * 1 . 5+(l  . 0-XX)**l  . 5)*CKLOG(l . 4142  14*SQRT  (XH)  + SQP.T 
1(1. 0-XX) )-l . 8856 188*XH**1 . 5-2. 828428*SQRT (XH)* ( 1 . 0-XX) ) *0.66666 
GO  TO  110 
END 


FUNCTION  E J 5 (3LZ  , BL  1 , BL2  , RRZ  , BR  1 , BR  2 , XI! ) 

C 

c 

c 

C EVALUATES  THE  EXPRESSION  FOP  THE  DEFINITE  INTEGRAL  OF  FTJX  OVER 

C TH  E LIMITS  XX-2H  TO  XX+2H  . BLZ  , PL  1 , EL  2 , ERZ  , BR  1 , 3R2  ARE  THE  COF.F- 
C FI  C I ENTS  OF  THE  QUADRATIC  APPROXIMATIONS  TO  P(X)  0 VER (XX-2H , XX ) AND 
C (X  X.XX+2H),  RESPECTIVELY. 

C 

C 

c 

C (5-PT.  INTKCRAL  ACROSS  LOG.  SINGULARITY.) 

C 
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FUNCTION  BINT3  (DUMl) 

REAL  KAY ( 50 ) 

COMMON  P(101 ) ,G(101 ) ,Q(101 ) ,H(101 ) ,X(10l ) ,XSQ(101 ) ,FX(101 ) ,ETNl  (10 
1 1 ) , ETN 2(101) ,ETP1 (101 ) , ETP2(101 ) , YYU1 (101 ) , YYU2 (101 ) ,TAU( 101 ) , ALP( 
2101) ,DEE(1G1) ,AA(50)  , KA Y , CA Y , DUM , A 1 , A 2 , AC A P , AL PH A , AMDA , AS Q , B B , B INT 
3 , CO, CL, CM, DELTA, DINT , ECAP, ELL , EM, EPS ILN , ETA  1 , ETA 2 , ETA  3 , TEE , ONE  , TJ  0 
4 , THR  , TIIET  A , TS  AVE  , XC  0 , XC  1 , XC  2 , XC  3 , XH  , XKZ  , NM  , NN  , JL  6 , JM6  , JR  6 , JL  1 3 , JM  1 
53, JR  13 , JTT 1 , JTM 1 , JTT  2 , JTM  2 , JTT  3 , JTM 3 , JTT  4 , JTM4 , JTT  5 , JTM  5 , JTT 7 , JTM  7 
6 , JTT 8, JTM 8, JTT 9 , JTM 9 , JTT 1 0 , JTM 1 0 , JTT 1 1 , JTM1 1 , JTT 12 , JTM 1 2 , JTT 21 , JTM 
7 2 1 , JTT 23 , JTM 23 , JTT 24 , JTM24 , JTT 25 , JTM 2 5 
NM-NN 

DO  100  1*1 , N M 

100  FX(I)*P(I )/ ( ELL -X ( I ) ) 

CALL  NT CRT E (BINT  3 ) 

RETURN 

END 


FUNCTION  BINT4  (DUMl) 

REAL  KA  Y ( 5 0 ) 

COMMON  r (101  ) ,G(101  ) ,Q( 101 ) ,H(101 ) ,X(101 ) ,XSQ(101 ) ,FX(101  ) , ETNl  (10 
1 1 ) , ETN  2(101 ) , ETP1 ( 10 l ) , ETP2 (101 ) , YYU 1 (101 ) , YYU2 (101 ) , TAU ( 1 01 ) , ALP ( 
2101)  ,DEE  (101)  ,AA(50)  , KAY  , CAY  , DUM  , A 1 , A 2 , ACAP  , ALPHA  , AMDA  , AS <) , BB  , BINT 

3,  CO, CL, CM, DELTA, DINT , ECAP, ELL, EM, EPS ILN, ETA  1 , ETA 2 , ETA 3, TEE, ONE,  TJ  0 

4 ,  THR .THETA , TS AVE , XCO , XC 1 , XC  2 , XC  3 , XH , XKZ , MM , NN , JL  6 , JM6 , JR  6 , JL 1 3 , JM 1 
53  , JR  13 , JTT  1 , JTM  1 , JTT 2 , JTM  2 , JTT 3 , JTM 3 , JTT 4 , JTM 4 , JTT  5 , JTM 5 , JTT  7 , JTM 7 
6, JTT 8, JTM 8, JTT 9, JTM9, JTT 10, JTM 10, JTT 1 1 , JTMJ 1 , JTT 1 2 , JTM 1 2 , JTT 2 1 , JTM 
72 1 , JTT 23 , JTM 2 3 , JTT 24 , JTM24 , JTT 2 5, JTM 2 5 

NM*NN 

DO  100  1*1, NM 

100  F X ( I ) = P ( I ) * S QRT ( X ( I ) ) /SORT (ELL -X( I ) ) 

CALL  NT  CRT  E (BINT4) 

RETURN 

END 
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FUNCTION  TEFUNC  (XNU , A , DELT  1 , DELT 2 , GAM  1 , GAM2 ) 
TFUN 


PLAYS  A PART  SIMILAR  TO  THAT  OF  ALFUNC. 


TEFUNC” (CAM 2 *XN  U+D  ELT  2 ) / (GAM  1 *XN U+D ELT l-SQRT (XN U * ( XN U + 1 . 0 ) ) ) 
TEFUNC =AT AN (T EFUNC ) -ATAN (DELT 2 /DELT 1 ) -ATAN (XNU /A) 

RETURN 

END 


FUNCTION  ALFUNC  ( XN U , A , OMGS Q , B ET 1 , B ET 2 , DE LT 1 , DE LT 2 , CAM  1 , CAM2 ) 
ALFN 


THIS  SPECIAL  FUNCTION  IS  PART  OF  EXPRESSION  FOR  F1(NU,A) — ALSO 
NEEDED  BY  F2(NU,A)  IN  THE  OFF-DESIGN  CALCULATION-. 


ALFUNC “ (CAM 1*XNU+DELT l-SQRT (XN U * (XN U + 1 . 0 ) ) ) * * 2+ (C AM  2 *XN U+D E LT 2 ) * * 2 
ALFUNC “AL  F UNC  * A*  * 2 / (DELT 1**  2+D  ELT  2**2)/(XNU**2 +A*  * 2 ) 

ALFUNC-CKLOG (ALFUNC) 

RETURN 

END 
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FUNCTION  BICM  (A) 

REAL  KAY (50 ) 

COMMON  P(101),C(101),Q(10l),R(101),X(10l),XS<)  (101), FX (101),  ETN  1(10 
li),ETN2(I01) , ET  P 1 (101) ,ETP2(i01) , YYU I ( 1 0 1 ) , YYU 2 ( 1 0 1 ) , TAU ( 1 0 1 ) , ALP ( 
2101 ) , DEE (10  I ) ,AA(50) , KAY , CAY , DUM, A1 , A2 , ACAP , ALPHA , AMDA , ASQ  , BB  , B INT 
3,  CO,  CL,  CM,  DELTA,  DINT  , EC  A P,  ELL  , EM , EPS  ILN  , ET  A 1 , ET  A 2 , ET  A 3 , TE  E , ON  E , T J 0 

4 , TIIR  .THETA  , TS  AVE  , XC  0 , XC  1 , XC  2 , XC  3 , XH  , XKZ  , NM  , NN  , JL  6 , JM6  , JR  6 , JL  1 3 , JM  1 

5 3 , JR  1 3 , JTT I , JTM 1 , JTT  2 , JTM2 , JTT 3 , JTM3 , JTT  4 , JTM4 , JTT 5 , JTM 5 , JTT 7 , JTM 7 
6 , JTT 8, JTM8, JTT 9 , JTM9 , JTT 1 0 , JTM 1 0 , JTT 1 1 , JTM1 1 , JTT 12 , JTM  12 , JTT 21 , JTM 
721, JTT 23, JTM 23, JTT 24, JTM24, JTT 25, JTM25 

NM-NN 

DO  100  1=1, NM 

100  F X ( I ) = P ( I ) *X ( I ) 

CALL  NTGRTE  (BICM) 

RETURN- 

END 


FUNCTION  BTX7  (XX) 

C BTX7 

C 
C 
C 

C THIS  SPECIAL  ROUTINE  EVALUATES  A DEFINITE  INTEGRAL  OVER  (0,1) 

C WHOSE  INTEGRAND  INVOLVES  P(X)  AND  A LOGARITHMIC  FUNCTION  OF  THE 

C VARIABLE  OF  INTEGRATION.  THIS  INTEGRAL  IS  NEEDED  BY  THE  ROUTINE 

C II II . 

C 

C 


REAL  KAY (50) 

COMMON  P(101)  ,C(101  ) ,Q(  101  ) ,11(101  ) ,X(101)  ,XSQ(  10  l ) , FX  ( 1 0 l ) , ETN  1 ( 1 0 
1 1 ) , ETN  2 ( 1 0 1 ) , ETP1 ( 101 ) , ET?2( 101 ) , YYU 1(101) , YYU2 (101 ) ,TAU (101 ) ,ALP( 
2 101 ) ,DEE (10 1 ) , AA(50) , KAY , CA Y , DUMA , A 1 , A 2 , AC A P , AL PH A , AMDA, ASQ, BB , BIN 
3T, CO, CL, CM, DELTA .DINT , ECAP , ELL , EM , EPS ILN , ETA  1 , ETA  2 , ETA  3 , TEE , ON E , T1 
40  , THR  , TH  ETA  , TSA V£  , XC  0 , XC  1 , XC  2 , XC  3 , XH  , XKZ  , NM , NN  , JL  6 , JM6  , JR  6 , JL  1 3 , JM 
513  , JR  13  , JTT  1 , JTM1  , JTT  2 , JTM 2 , JTT  3 , JTM 3,  JTT  4 , JTM4  , JTT  5,  JTM.  5 , JTT  7 , JT ft 
67 , JTT 8, JTM 8, JTT 9, JT M9 , JTT 1 0 , JTM 10 , JTT 1 1 , JTM1 1 , JTT  12 , JTM  12 , JTT  21  , JT 
7 M 2 1 , JTT23, JTM2  3, JTT 24, JTM24, JTT2  5, JTM 2 5 
COMMON  /BLK1 /VlBAR, VI BARI , ECAP1 , YC  12 


C 


I 

-66-  16  June  1977 

JF: jep 

HYDROFOIL  DESIGN  PROGRAM 


G7(S)=ALOC(SQRT (ELL-XX) *SQRT (S)+SQRT (XX ) * S QRT ( ELL-S ) ) 
XCO-O. 0 
XC 1 *X ( NN -2 ) 

XC2«X(NN-1  ) 

XC3=X(NN) 

ETA 1 =G  7 (XC  1 ) 

ETA2=G7 (XC2) 

ETA3=G  7 (XC  3) 

CALL  PARAB  ( B3 0 , BB 1 , BB  2 ) 

NM-NN-2 
FX ( 1 ) = 0 . 0 
DO  100  1=2,  N'l 

100  FX(I)=P(I)*G7(X(I)) 

CALL  NTGRTE  (TGI) 

NM  = 3 
J-NN-2 

DO  110  1=1, NM 

FX(I)=C(J)*C7 (X(J)  ) 

J=J  + 1 

110  CONTINUE 

CALL  NTGRTE  (TG  2 ) 

BTX7=TG 1+TG  2+BB*F0UR (X(NN-2 ),BB0,BB1,BB2) 

RETURN 

END 


Sample  Calculations  for  Estimation  Mode 


16  June  1977 
JF:jep 


© 

• r* 

O' 

o o 

n* 

r* 

• 

it 

o 

■ *r 

«4  • 

BS  O 

• 

A« 

*9 

«3 

u 

4 ■ 

40 

• ■ 

O k> 

M 

M O 

O U 

O 

o 

O M 

O N 

O 

O ft* 

M • 

O M 

r*  fti 

O iJ 

<*)  VO 

o 

n n 

in 

o m 

• *• 

m cd 

o 

o o 

OS 

• 

m r» 

ft] 

ft] 

• r» 

00 

• o 

o • 

O' 

■c  o 

o 

r* 

o o 

O' 

B O 

■ 

«<  o 

es  N 

m 

o 

«<  CD 

04 

O O 

© 

• 

as  m 

X 

o 

© • 

o o 

« 

N 

r* 

as  N 

o 

O — 

ni  i 

*-4 

M 

o w 

ft] 

H -< 

i <n 

O U 

ft]  r* 

ft* 

as  wo 

m oo 

o 

© o 

o CO 

ft  «*- 

in  o' 

1- 

• 

04  O 

os 

• O 

CO  9 

o 

o 

•“  • 

<0 

as 

o o 

IN  ft 

• 

o as  o 

o 

1 o x 

9 

ta  m 

a: 

o H 

ii  m 

04 

o «< 

Q 

o 

o h 

O L> 

1 

O M 

O 

ft] 

O > O 

o o 

NO 

O < o 

»•  O 

O U *“ 

* O *“ 

03 

o • 

O r-  o 

in 

• O 

• **  1 

NO 

o 

n o o w 

NO 

ioooooooooooooooooo 

• I « I I I I I I I I I I I I • I 
IUUUUUUUUUWUUUUUUUU 
ir^cn<Nco^nn*-orj»“oor'(N^cD^ 
lOcoiriyflr^cocDvona'iOvflrnonna' 
,rH0'r>u^vor^^r^^'iAnr9n^)^r"C'r'j 
i co  n»  *—  vo  O'  co  ojo'nipcrcT'Qujino' 

>rna^-r*^r-i^r»o^o*ocnocT'r'sou-> 
ioa»voNO'»onocoinnocDofn*-o'r' 
• *-a'cr»(^co(Dco=on-r^r^r^vow3vov^u*>m 


OOOOOOOOOOOOOOOOOOOOOOOOOOOOOO 


OOOOOOOOOOOOO' 
1111(1111111 
UDJUUUKMHKUUUU! 
* O3Df**r'»C0*“CNUDOr^m<Nv£!< 
oc'C'Nr*tnu^rN^fNff'^r-. 
OU’>r,*OvDvOvDl/>vOOOOvO< 
Q«c'v“r^cyiT^crcTr^ir»C'f-c'< 
t^c'C'^rvnma-r^vDvDCT'cootN. 

wo'C'3hnoNin»-c'(7'»-n| 

C'^CM»-»-»-CCvDir.  rn<NCS|r-. 


I M (M  N (M  ( 

» o o o o c 

(III 

I w u U W I 

i 04  O'  CO  Ol  v 

i cn  r-  m o c 

i lO  04  <N  O'  : 

*-  m rg  cc  c 
10505^1 
04  04  04  4*3  r 


I 04  04  (M  04  < 

* O O O O 1 

I I I I 

I U u U (O  i 

• cy  vo  st  o < 

> m *-  04  < 

I r n pi  O I 

■040®  Oi 

) o ® y i 

) y o o i 

i y y y i/i  i 


O 1* 

x ||  ft  y • 

Ii  o cn 

in  o'  co  o 

S3 

NO  it 

X ft] 

04  nO 

o 

o«n  ••  h 

ft]  o 

in  st  in  *- 

N O 

o m . z ft 

MOO 

O O 04  o O •< 

in 

in  o o 

O * • M U 

at 

• O 

WOO  H w 

o 

ft  o ft 

vo  «-  i n «c 

ij 

ft]  o o 

M • 0-  *J 

ij 

f-»  • o 

04  O If  •<  O 

o 

in  o o 

r-  C U ’J 

ft. 

o 

®C  41  u < w 

o 

O f'  H ^ »- 

-< 

D O 

1 »-  u o 

h 

• o 

ft)  4 1 

0t 

a.  *- 

o o as  co 

o 

in  c # 

r*  o oo 

o o 

O'  CM  *-  M rr\ 

H 

H U 

04  ii  r-  o o in 

© 

as 

IN  D O I CO  04 

ft 

* 1-  M 

CO  SC  4 CO  w O O 

as 

o 

y < w y y » *- 

M 

• o w •- 

r-  CC  O'  ft.  *- 

cn  os  to  o 

• O'  O'  m e».  • 

ft. 

h •- 

o o*  ry  r»  o o 

O 

tft  • ^ W 

O'  O'  o 

1 O o 

O'  ft  n a 

O 

© in  m o' 

II  O O'  r»  r-  o M 

a*. 

M O' 

O 0«  v • ft.  cn 

M 

es  ft*  « o 

<<  in  • o o < 

DOW® 

WOO  H CD 

in 

t-  r» 

t © *“ 

M 

ft  • oj  n 

a:  o n it  o.> 

►1 

O ft-  n 04 

m h o •-  a: 

ft]  •<  • 

n •-)  \ u m 

• 

• O ft:  o 

os  as  u M x 

O U -a 

o 

as  ft*  It 

fti 

•< 

a 

ft] 

o o o o 

I I I 


I O'  04  04  OI 

ft)  O O O O 

O I I I I 

st  U U U U 

CO  O 04  04  lO 

m c n n o 

vo  rn  r*  m co 

vo  <mvozro* 

m H r>  m vo 

l W\C  CC  3 r 

o cr  •-  oi  cm 


0 o © o o « 

1 I I I I 

04  04  04  04  04  < 
O O O O O « 
I « I I ( 
u u u u U I 

y o®  *- 

»-  ® y oio  « 
ffl  *-y  o<o« 

DOifltf  n 
o •*  o y vc  i 
*-  cr  m C"  co  < 
n n n n n • 


000000004 

• I • I I • 

01040404040404041 
OOOOOOOO. 
I I | • » I I I 

MUUUUUUUI 

mono^  voop»i 
C'fv0p*00l004. 
nooCUT®  r-vfl« 
040404\0'0040'0 
lT040004n04f-< 

OOO^^Off-i 


•oooooooooooo 


I 04  04  04  04  ( 

> O O O O < 

I I # » 

IUUUUI 

'O  no  a i 

IO**f  ji 

> CD  vO  VO  04  l 
•O  rn  O'  O'  I 
l O'  M y vO  l 
1 CO  C3  O O l 
in  n n n 1 


4 04  04  04 

> O O O 

I I I 
IUUU 
l 04  C 04 
N nD  O O 

> O'  O O 

> 04  O 04 
' T»  CO  CO 

-noj»- 
i n n n 


©OOOOOOOO 
I I I I I I I I I 


OOOOOOOOOOOOOOOOOOOOO 
I I I I I I I I I I I I I I I I I I I I I 


®0404*-»-*-^ff- 
OOOOOOOOO 
I I I I I I I I I 
u^uwwwwuw 
o in  o y rr  o*  co  04  •- 
04  04*-r-unr00in04 
in  no -o  n- f-fv 

CDOONOr-f^nO 
►-inp^oo^r-comr-o 
vDonc4inooo4in 
(N  in  O *-  r-  r-  04  N 04 


OOOOOOOOOOOOOOOOOOOOO 
I I I I I I I I I I I I I I I I I I I I I 
UUUUUMMMblUUMUUUtUUUUMU 
infflinc'ooinoo4®<70sino«f'innc4ffl 
fsyNyoinnnvonfNO'inyON'-oiNOO 
®n^mo'*oon(M®(NOn»-np>o'offlC 
04na'Offl04nfioa'C'03piiN'nNohfNo 
nmoffiffio^onMomycN^OMnoioo 
ooi*-nL''oc'^*f,>inoc'^ninocoojyin 
N04nnnr*>n^yy  yyinminininoooo 


OOOOOOOOOOOOOOOOOOOOOOOOOOOOOO 


> O O o 

I I I 

) ft]  ft)  to) 

> o o CO 

> O o o 

> O o O' 

> O O O' 

> O O O' 

> o o O' 
4 y in  n 


^-OOOOOOOOOOOOOOOOOOOOOOOO 

OOOOOOOOOOOOOOOOOOOOOOOOO 

I 

UUMUUUUUMUUUUCUUUMUUUUUCUUU 

VDC'OOO'OOO'OOOOOOOOOOOOOOOO'O 

O'O'OOO'OOO'OOOOOOOOOOOOOOOC'O 

O'O'OOO'OOO'OOOOOOOOOOOOOOOC'O 

O'C'OOC'OOO'OOOOOOOOOOOOOOOO'O 

O'C'OOCT'OOO'OOOOOOOOOCOOOOOO'O 

C''*y®noNniflK)on4yvooOn4jJo®oo4yin® 

cr^ffr-fNMiNNiNnnnnnyyyyyminininin 


OOOOOOOO©' 


i o o o O o I 


lflr’ifl,“'Or,'0'“'Or“i4)^,®»"i4>»*NOfiOr“'0 


Sample  Calculations  for  Design  Mode  (Cont. 


T 


-70- 


16  June  1977 
JF:jep 


0 o o o o 

1 I I I I 
M M W W W 

r-  n a in  vO 
vo  n o'  o in 
i voo' 
in  in  m in  m 
inn*-o'h 
in  m in  ^ ^ 


0 o 

1 i 
w co 

o 

VO  vO 

o> 
cm  cr 
X.  VO 
in  m 
a st 


0 o o o o 

1 I I I I 
co  co  w w w 

NfMP^Nr- 

c>  ^ O'  m »• 
n o o'  n O' 
sOvONvO  ^ 

o o vo  in  ^ 

r-  o'  r»  in  n 

?nnnn 


o o o o 
till 

uuuu 
co  in  oo  cm 
vo  o st  cn 
r-  PM  ^ *- 

n*  ^ n 

CN  O vO  r— 
»•  O'  lO 
CO  CM  CM  CM 


0 o o o o 

1 I I I 

uuuuu 

vO  tfv  O O 
co  o ctv  r*  o 
co  o o 
•-  o o 

^ OVOO 

r*CO  o o o 

N«-  r-f-O 


o 

o 

O 

o 

O 

o 

o 

O 

O 

O 

o 

O 

o 

© 

O 

o 

o 

© 

o 

o 

o 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

o 

O 

O 

o 

o 

O 

O 

O 

O 

O 

O 

O 

o 

O 

O 

O 

O 

o 

O 

O 

o 

1 

1 

1 

1 

1 

1 

1 

» 

1 

1 

8 

1 

1 

1 

• 

1 

1 

1 

I 

1 

• 

co 

w 

co 

co 

w 

co 

CO 

CO 

CO 

CO 

W 

CO 

CO 

CO 

CO 

CO 

CO 

w 

CO 

CO 

CO 

vO 

*” 

CD 

o 

vO 

O' 

o 

vO 

CD 

O' 

in 

m 

CM 

CO 

CM 

CO 

o 

9 

CO 

•— 

co 

vO 

IT. 

o- 

O' 

CM 

VO 

vO 

PM 

<*> 

CO 

vC 

ST 

VO 

CM 

CO 

CM 

vO 

o 

o 

o 

cr 

O' 

O' 

CD 

CM 

CM 

co 

CO 

O' 

O' 

O 

in 

tn 

co 

O' 

n* 

sO 

O' 

r— 

o 

r- 

cr 

o 

o 

ST 

m 

VO 

■o 

vo 

c- 

CO 

o 

ST 

in 

CM 

r» 

O' 

CO 

VO 

O' 

*— 

CM 

co 

o 

CM 

f— 

O' 

r- 

+— 

00 

o 

in 

o 

in 

O' 

CO 

r- 

r- 

O' 

CM 

ST 

vO 

CO 

O 

CM 

CO 

in 

O' 

o 

CM 

m 

r* 

CO 

O' 

r— 

CM 

m 

m 

VO 

vO 

vo 

vo 

r* 

P- 

r- 

p- 

p- 

p- 

CO 

CO 

00 

CO 

CD 

CO 

CO 

O' 

O' 

o 

o 

O 

O 

o 

o 

o 

o 

o 

o 

o 

o 

O 

o 

o 

o 

O 

o 

o 

o 

o 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

i 

CO 

CO 

m 

CO 

co 

VO 

O 

o 

O 

O 

O 

O 

O 

O 

o 

O 

O 

O 

O 

O 

o 

o 

o 

o 

o 

o 

o 

| 

1 

• 

1 

» 

1 

1 

1 

1 

1 

1 

t 

1 

1 

1 

I 

I 

i 

1 

1 

M 

CO 

co 

« 

CO 

co 

CO 

CO 

CO 

CO 

CO 

CO 

CO 

CO 

CO 

« 

CO 

pj 

CO 

CO 

CO 

st 

CM 

co 

O' 

O' 

O' 

r- 

in 

CO 

CM 

vO 

ST 

vO 

CM 

in 

ST 

CO 

o 

o 

St 

O 

O' 

o 

m 

r— 

in 

in 

co 

r« 

vO 

O 

r* 

CO 

r* 

PM 

CO 

CO 

r— 

VO 

m 

CM 

vC 

O 

ST 

o 

O' 

CM 

O 

CO 

c-» 

r- 

VO 

o 

o 

CM 

o 

CO 

r- 

O' 

e- 

O' 

ro 

CM 

CO 

O' 

r- 

o 

CO 

r- 

*— 

CO 

r- 

V- 

in 

O' 

9 

o 

CO 

vO 

in 

CM 

o 

r* 

n 

O' 

in 

r— 

vo 

o 

O' 

CM 

VO 

CM 

O' 

CO 

in 

r* 

O 

O' 

CD 

p* 

in 

CM 

r— 

O 

CO 

p^ 

in 

CO 

CM 

o 

O' 

*— 

ST 

vO 

CO 

VO 

r> 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

CM 

r- 

r» 

c- 

CO 

in 

co 

r~ 

r* 

ooooooo  ooooooooo ooooo 

• I ••  I I f I I I I I I I • I 9 111 


1 


DISTRIBUTION  LIST  FOR  UNCLASSIFIED  TM  77-199  by  J.  Fernandez,  dated 


16  July  1977 

Commander 

Naval  Sea  Systems  Command 
Department  of  the  Navy 
Washington,  DC  20362 
Attn:  Library 

Code  NSEA-09G32 
(Copy  Nos.  1 and  2) 

Naval  Sea  Systems  Command 
Attn:  A.  Bress 

Code  NSEA-03 
(Copy  No.  3) 

Naval  Sea  Systems  Command 
Attn:  L.  Benen 

Code  NSEA-0322 
(Copy  No.  4) 

Naval  Sea  Systems  Command 
Attn:  T.  E.  Peirce 

Code  NSEA-0351 
(Copy  No.  5) 

Naval  Sea  Systems  Command 
Attn:  D.  Creed 

Code  NSEA-03132A 
(Copy  No.  6) 

Commander 

Naval  Ship  Engineering  Center 
Washington,  DC  20360 
Attn:  R.  J.  Cauley 

Code  NSEC-6140B 
(Copy  No.  7) 

Naval  Ship  Engineering  Center 
Attn:  F.  Welling 

Code  NSEC-6144 
(Copy  No.  8) 

Commanding  Officer 
Naval  Underwater  Systems  Center 
Newport,  RI  02840 
Attn:  Library 

Code  LA15 
(Copy  No.  9) 

Naval  Underwater  Systems  Center 
Attn:  J.  D.  Powers 
Code  SB34 
(Copy  No.  10) 
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Commanding  Officer 
Naval  Ocean  Systems  Center 
San  Diego,  CA  92132 
Attn:  A.  G.  Fabula 
Code  5002 
(Copy  No.  11) 

Commanding  Officer  & Director 
David  W.  Taylor  Naval  Ship  R&D  Center 
Department  of  the  Navy 
Bethesda,  MD  20084 
Attn:  W.  E.  Cummins 
Code  15 
(Copy  No.  12) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  S.  F.  Crump 
Code  1505 

(Copy  Nos.  13  - 32) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  R.  Wermter 
Code  152 
(Copy  No.  33) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  W.  B.  Morgan 
Code  154 
(Copy  No.  34) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  R.  Cumming 
Code  1544 
(Copy  No.  35) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  J.  McCarthy 
Code  1552 
(Copy  No.  36) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  T.  Brockett 
Code  1544 
(Copy  No.  37) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  Y.  T.  Shen 
Code  1524 
(Copy  No.  38) 
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David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  M.  Sevik 
Code  19 
(Copy  No.  39) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  W.  Blake 
Code  1942 
(Copy  No.  40) 

David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  T.  C.  Mathews 
Code  1942 
(Copy  No.  41) 


David  W.  Taylor  Naval  Ship  R&D  Center 
Attn:  Tech.  Info.  Lib. 

Code  522.1 
(Copy  No.  42) 


Commanding  Officer  & Director 
David  W.  Taylor  Naval  Ship  R&D  Center 
Department  of  the  Navy 
Annapolis  Laboratory 
Annapolis,  MD  21402 
Attn:  J.  G.  Strieker 
Code  2721 
(Copy  No.  43) 

Commander 

Naval  Surface  Weapon  Center 
Silver  Spring,  MD  20910 
Attn:  J.  L.  Baldwin 
Code  WA-42 
(Copy  No.  44) 

Office  of  Naval  Research 
Department  of  the  Navy 
800  N.  Quincy  Street 
Arlington,  VA  22217 
(Copy  No.  45) 

California  Institute  of  Technology 
Division  of  Eng.  and  Appl.  Sci. 
Pasadena,  CA  91125 
Attn:  A.  J.  Acosta 
(Copy  No.  46) 


California  Institute  of  Technology 
Attn:  T.  Y.  Wu 
(Copy  No.  47) 

Dr.  B.  R.  Parkin 

The  Pennsylvania  State  University 
APPLIED  RESEARCH  LABORATORY 
Post  Office  Box  30 
State  College,  PA  16801 
(Copy  No.  48) 

Mr.  J.  Fernandez 

The  Pennsylvania  State  University 
APPLIED  RESEARCH  LABORATORY 
Post  Office  Box  30 
State  College,  PA  16801 
(Copy  No.  49) 

GTWT  Library 

The  Pennsylvania  State  University 
APPLIED  RESEARCH  LABORATORY 
Post  Office  Box  30 
State  College,  PA  16801 
(Copy  No.  50) 

Defense  Documentation  Center 
5010  Duke  Street 
Cameron  Station 
Alexandria,  VA  22314 
(Copy  Nos.  51  - 62) 
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